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The problem of a transformer of 7 windings, each con- 
nected to a circuit containing a source of electromotive 
force, is treated. The assumptions are pure reactance 
windings, perfect couplings and linear loads. The equations 
for the current in any circuit, the open-circuit voltage 
across any winding, the short-circuit current in any 


INTRODUCTION 
N linear electric circuits containing trans- 
formers of only a few windings, the various 
currents and voltages are readily obtained by 
applying Kirchhoff’s laws to the given circuit or 
to some equivalent circuit. However, when 
circuits containing transformers of a large number 
of windings are considered, this method becomes 
very laborious; the amount of labor being 
proportional to the factorial of the number of 
windings. Since the Kirchhoff law method is, in 
general, the only way in which circuits of this 
type can be solved, it is unfortunate that the 
necessary work increases at such a rapid rate; 
nevertheless, there seems to be no way around 
this difficulty except by making certain simpli- 
fying assumptions relative to the transformers. 
There are two assumptions which are usually 
made in treating transformers. In the first place, 
it is assumed that the mutual impedance between 
any pair of windings is equal to the square-root 
of the product of the self-impedances of the 
windings;' secondly, it is assumed that the self- 
impedances of the windings are very large in 
‘ This involves assumptions of pure reactance windings 
and perfect couplings. Timbie, Elements of Eelectricity, John 
Wiley and Sons, pp. 380-383 (1925). 
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winding, the impedance looking in at any winding, the 
voltage transformation ratios, the current transformation 
ratios and the impedance transformation ratios are given. 
The ideal n-winding transformer, which is a more re- 
stricted case, is also treated. 


comparison to the load impedances. A trans- 
former satisfying both of these conditions is 
known as an ideal transformer.” 

In most actual transformers, so far as their 
transmission properties are concerned, very little 
error is introduced by these assumptions. They 
are both involved in some of the best known and 
most used two-winding transformer relations, 
such as the ratio of the primary current to the 
secondary current being equal to the ratio of the 
secondary turns to the primary turns. The latter 
is probably the more damaging of the two 
assumptions since transformer circuits are fre- 
quently found in which the load impedances are 
high in comparison to the impedances of the 
transformer windings. This is exemplified by the 
fact that the loads on input transformers in 
vacuum tube circuits are often extremely high. 

By using only the first of these assumptions, 
sufficient symmetries are introduced into the 
determinant of the Kirchhoff law equations to 
enable the currents, voltages, etc., to be found for 
even the general case of a transformer of 
windings. The formulae for the ideal n-winding 

*K. S. Johnson, Transmission Circuits for Telephonic 


Communication, Chapters VI and VII, D. Van Nostrand 
(1927). 
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transformer are obtained from these by making 
use of the second assumption. 

The mathematics necessary for solving such a 
set of m simultaneous equations is somewhat long 
and involved and has a tendency to camouflage 
the results so it has been omitted. The author 
will be glad to furnish any one writing him a 
complete proof of all of the formulae given here. 


FORMULAE 


Consider the n-winding transformer and con- 
nected circuits shown in Fig. 1.4 The terminals of 
the transformer are so numbered that all 
currents flowing from odd to even numbered 
terminals give magnetomotive forces in the same 
direction in the core. The impedances are, in 
general, complex quantities. The electromotive 
forces are all of the form | E|e” where @ is the 
angle of E in the complex plane using, say, the 
maximum electromotive force from the odd to 
the even numbered terminal of the corresponding 
circuit as a reference. 
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Fic. 1. Circuit diagram. 


The current in the pth winding, where p is any 
number between one and 2, is 


E, Zw VEZ 3 /Z,! 
[,= oe. 2 a (1) 


4) ZNA4FD Z/Z,’) 
k=I1 


* The circuits are shown as consisting of simply a source 
of electromotive force in series with an impedance; this 
is perfectly general since, according to Thevenin’s theorem, 
any two-terminal circuit can be reduced to this form. See 
M. L. Thevenin, Sur un nouveau theoreme d'electricite 
dynamique, Comptes Rendus 97, 159 (1883). Also K. S. 


Johnson, reference 2, p. 79. 
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The open-circuit voltage in the pth circuit 


EF, +> Zi Z;')—Z,,’ > a EZ) /Z). 
(Vo) p= k=1 = kel (2) 


1+>- Z; Z:-Z, > Py 
k=1 





The short-circuit current in the pth circuit is 


(1.),)=(E,/Z,)A+E Z,/Z,') 
k=1 


h 


nm 
“An *Dd FiZit/Z,'. (3) 
A l 


It may be noted that if all the electromotive 
forces except one, for example E,, are zero, the 
short-circuit current in the pth circuit is 


(1,) p= —(Eq/Zq)\(Z_/Zp)'. (4) 


This current depends only on the constants of the 
pth and gth circuits. 

The impedance looking in at the pth winding 
(i.e., the impedance as measured across the two 
terminals of the pth winding) is 


(Z)p=Zp/(1+D. Ze/Ze'—Z,/Z,"]. (5) 
k=1 


The voltage transformation ratio (i.e., the 
ratio of the voltage applied across Z, to the 
voltage across Z,’ due to this voltage across Z,,) is 


(V,/V,)=—(Z,/Z,). (6) 


The current transformation ratio (i.e., the 
ratio of the current in the pth circuit to the 
current in the gth circuit due to the current in the 


pth) is 


(1/1) = —Z (ZZ) TAFE Ze/Ze! 
r=] 


A 


—Z,/Z,')\. (7) 


The impedance transformation ratio (i.e., the 
ratio of the impedance looking in at the pth 
winding to the load impedance in the gth 
circuit) may be found by dividing Eq. (5) 
through by Z,,’. 
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IDEAL’ TRANSFORMER FORMULAE 
The current in the pth winding is 
n 
E Z,* > E,Z). ta 
. k=1 
I, = - —e (3) 


A Ai EA 


k=1 





The open-circuit voltage in the pth circuit is 
n n 
E, >, Z:/Z:'—-Z,' > ExdZ,'/Z,' 


, k=1 k= 
(v),.——— : ~~ 


a Ae ee, 


A=1 





The short-circuit current in the pth circuit is 


(1,),=(E,/Z >) >. Z, Ze 


k=l 


—Z,°> E,Z,3 Z,{. (10) 
k=1 


The impedance looking in at the pth winding‘ 


is 
(Z),=@Zp/l >, Zt/2Z: —Z,/2,' |. (11) 
The voltage transformation ratio is 


(V,,/V.) = —(Zp/Z,)}. (12) 


The current transformation ratio is 


(1 p/Lq) = —Za (ZZ 4)" > 2/24’ —Z)/Z,'). (13) 


~ 


The impedance transformation ratio may be 
found by dividing Eq. (11) through by Z,’. 


4 This equation, for p=1, was arrived at inductively by 
K. S. Johnson and is found, in different form, in his book, 
reference 2, p. 41. The symmetry of this relation suggested 
the possibility of treating the m-winding transformer. 
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Methods for Calculating the Volumetric Composition of Fluid Mixtures* 


S. H. INGBERG, Bureau of Standards 
(Received November 9, 1933) 


The mathematical treatment relating to mixture of fluids 
in a container maintained at constant pressure and volume 
is developed, assuming perfect and instantaneous diffusion, 
free efflux and no reaction between the components of the 
mixture. For comparison, the equations for mixture within 
a container from which no efflux takes place are given and 
also a treatment applicable to intermediate conditions. 


I. INTRODUCTION 
WO limiting conditions under which fluids 
are mixed are here discussed, one for which 
the volume and pressure of the fluid is maintained 
at initial values by providing means for free 
efflux, and the other involving retention of all 
original and added quantities during the mixing 
operation. The latter condition is possibly more 
often present and the methods of computing the 
volume of components are more obvious. The 
methods corresponding to the former condition 
do not appear to be extant and those applicable 
for the latter have been used instead, sometimes 
with considerable error in results. 
The mathematical treatment developed applies 
in general to mixtures of fluids in a container or ° 
inclosure maintained at constant pressure and 
volume, assuming perfect and instantaneous 
diffusion, free efflux, and no reaction within the 
mixture. For comparison, the equations for 
mixture within a container from which no efflux 
takes place are given and also a treatment 
applicable to intermediate conditions. Specific 
applications are developed pertaining to the use 
of inert gas, such as carbon dioxide, released 
either from the liquefied form or as one of the 
inert components of flue gas or engine exhaust 
gas, to prevent fire and explosion in connection 
with hazardous processing and storage. Either an 
inert atmosphere is maintained continuously 
within the enclosure or a quantity of gas is stored 
* Publication Approved by the Director of the Bureau of 

Standards of the U. S. Department of Commerce. 


Specific applications are developed pertaining to conditions 
involving inert gas, such as carbon dioxide, flue gas or 
engine exhaust gas, released into an inclosure by automatic 
or manual means to extinguish fire, and to the subsequent 
ventilation of spaces thus deluged to reduce the toxic gas 
content and raise the oxygen content to limits making such 
spaces safe for entry. 


ready for release by manual or automatic means 
in case fire occurs or a condition develops 
conducive to fire or explosion. The present 
mathematical treatment applies for the latter 
condition and to requirements for subsequent 
ventilation of the spaces thus deluged to reduce 
the toxic gas content and raise the oxygen 
content to limits making such spaces safe for 
entry. 


II. ADDITIONS OF A PURE FLUID REQUIRED FOR 
A DESIRED CONCENTRATION 


1. Free efflux 


If V the volume or quantity of a fluid initially 
present in a container and y volume or quantity 
units of another fluid, each equal to V, are 
admitted to give a concentration x of the fluid 
added, then the change in concentration dx 
produced by the admission of a small increment 
Vdy of the added fluid will be proportional to the 
difference between this increment and _ that 
contained in the mixture forced out of the 
inclosure by the admission of Vdy. Accordingly 
assuming an ideal condition involving perfect and 
instantaneous diffusion and also 
within the mixture, 


no reaction 


dx =(Vdy—xVdy)/ V 
or dy=dx /(1—x); hence 
y= —log, (1-—x)= 


— 2.303 logio (1 —x).! (1) 


! A previous derivation of this formula by E. R. Weaver, 
Report No. 40, National Advisory for Aeronautics, is noted. 
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CONCENTRATION, X 


Fic. 1. Relation of volumes of fluid admitted to its con- 
centration 


The plot, Fig. 1 (lower curve) shows the 
general relation between the number y of 
volumes of fluid that must be admitted to attain 
a concentration of x within the receiving 
container, assuming the rate of efflux to equal 
that of the admission. 

If an initial concentration x; of the fluid being 
added is present, the requirements for increasing 
the concentration to x can be obtained from the 


relation, 
1-—x 
a 
1 om 3s 


= — 2.303 log io 


dx 


. -/ te 
1-x 

| (2) 
1-—x, 


The conditions for which the above and related 
subsequent developments hold are fulfilled where 
a gas or liquid is admitted into a container, 
closed except for inlet and outlet openings, that is 
filled respectively with a gas or liquid mixture, 
with conditions for diffusion and efflux con- 
forming with the assumptions made, or where 
liquids are mixed under similar conditions except 
that the container is open with means for efflux 
such that the volume of the mixture will not 
exceed the original volume present. 
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2. No efflux 


If no efflux takes place while the fluid is 
admitted, the concentration x after y volume or 
quantity units are admitted will be 


x=y/(1+y), 


The upper curve in Fig. 1 is a plot of this 
equation. By comparison it is seen that for a 
given concentration a greater addition is required 
than for the condition of free efflux. 

For an initial concentration of x; 


x=[(xit+y)/(i+y)] y=((x—*1)/(1—x) ]. 


This gives the requirements for increasing a 
concentration to x from an initial value of x. 


or y=x/(1— x). 


or 


3. Intermediate conditions of efflux 


While the conditions approximating free efflux 
and no efflux are generally more readily estab- 
lished where mixture of fluids is involved, a 
method whereby determinations can be made for 
intermediate conditions may be of interest. The 
one intermediate condition for which the 
derivation is made, involving a constant ratio of 
efflux to influx, will not in general obtain, but a 
method is indicated whereby it can be applied 
where this relation varies. 

Assuming a ratio R of rates of efflux to influx 
into a chamber containing initially a quantity V 
of a fluid, the addition of an increment Vdy of a 
component x of the mixture will change its 
concentration according to the relation, 


dx = Vdy/| V[1+y(1—R) }} 
~—xVdy/{V[1+y(1—R)]}. 


The first term of the second member represents 
the portion of the increment Vdy that is mixed 
with the original amount V of the fluid, and the 
second term, the part of the same increment that 
along with other components is forced out of this 
original volume or quantity of fluid by the 
admission of the given portion of Vdy. Reducing, 


dy/(1+y(1—R) ]=dx/(1—x) 
whence 
log.[1—y(1—R)]/(1—R) = —log. (i—x). (3) 


This can also be written as, 


log.[1—y(1 — R) ]/(1—R) = —log, (1—x), 
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since the equality will hold for logarithms taken 
to any base. For free efflux, R=1 and Eq. (3) 
becomes 


y = —log. (1—x) (1) 


for R=0, Eq. (3) reduces to y=x/(1—x), which 
expresses the relation for no efflux. 

The intermediate curves in Fig. 1 are plots of 
Eq. (3) for four values of R. 

Where R varies, successive computations can 
be made of the concentration from suitably 
spaced simultaneous readings or determinations 
of influx and efflux, assuming R for the interval 
between readings to be constant. For each 
determination, the increment y between intervals 
is computed as the fractional part of the quantity 
of fluid present at the beginning of the interval 
and the computation will be for the total 
concentration at the end of the interval. For this 
purpose the differential equation from which Eq. 
(3) is derived is integrated between the concen- 
tration limit x, at the beginning and x at the end 
of the interval to give, 


loga[1+¥(1—R) ]/(1—R) 
= —log, [(1—x)/(1—-1) ]. (4) 


With suitable spacing of readings for influx 
and efflux, the error due to the assumption of a 
constant ratio between them for each interval 
can be made so small that it will not affect the 
order of accuracy obtainable for the given type 
of determination. 

Il]. AppITIONS REQUIRED OF A MIXTURE OF 
FLuips TO GIVE A DEsIRED CONCEN- 
TRATION OF ONE oF ITS 
COMPONENTS 


1. Free efflux 


If the fluid being added contains the fractional 
part x2 of the component whose concentration is 
to be changed, for the addition of a small 
increment Vdy of this fluid mixture under 
conditions of free efflux, the change dx in the 
concentration of the given component will be, 


dx = (x2Vdy—x Vdy)/ V =(xe—x)dy, 
whence, 


y= —log, (x2—x) = —2.303 login (xe—x). (5) 
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For an initial concentration x; of the given 
component in the mixture to which the addition 
is made, the following relation is obtained: 


* dx Xe—-X 
y= { —— = —log, |= 
V2, Igo X xe X1 


Xe—~X 
= — 2.303 logiy || (6) 


te— 23 
This can also be written in the form, 


xX =X2—(xe—x1)e’. (6’) 

The above can be taken as expressing a general 
relation from which the special cases for free 
efflux previously discussed can be derived. Thus 
for the addition of a pure fluid, x. =1 and Eq. (6) 
reduces to Eqs. (1) or (2), depending on whether 
or not an initial concentration of the component 
being added is present in the fluid to which the 
addition is made. 

Whether an increase or decrease in the 
concentration of a component is to occur will 
depend on the magnitude of its concentration x2 
in the fluid being added, relative to its concen- 
tration x or x; in the fluid to which the addition is 
made. In general as the addition continues x will 
approach x2 in value. For unlimited additions of a 
pure fluid, that to which the addition is made will 
become identical with it in composition. At any 
stage if xe is greater than x, increase in concen- 
tration will occur, if smaller, dilution. For the 
case x2=0, none of the component concerned is 
present in the fluid being added and with 
continued additions its concentration in the 
fluid to which addition is made will approach 
zero. 


2. No efflux 


If a concentration x2 of the component whose 
concentration is to be changed is present in the 
fluid added, the concentration x in the fluid to 
which y volume or quantity units is added, each 
equal to that originally present, is given by, 


x=xey/(1l+y) or y=x/(x2e—x). 


If an initial concentration x, of the component 
concerned is present in the fluid to which the 
addition is made, 
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x= (xi txey)/(1+y) 


The above can also be regarded as a general 
equation from which the relations for special 
cases of increase of concentration or of dilution 
can be derived for the condition of no efflux. 


or y=(x—4%x)/(xe—x). 


3. Intermediate conditions of efflux 
For the condition that the fluid admitted 
contains the fractional part x2 of the component 


whose concentration is to be increased, Eq. (3) is 
modified to give, 


log. 1+9(1 — R) ]/(1—R) = —loga (v2—x). (7) 


Where an initial concentration x, of the 
component is present in the fluid to which the 
addition is made, Eq. (7) changes to, 


log.l1+(1—R)]/(1—R) 


= —log. [(x2—x)/(xe—4%1) ]. (8) 
Eq. (8) can also be written in the form, 
y= t1/(1—R)} {L@2—41)/(x2—-x) PF -1}.  (8’) 


The above equations express in their most 
general form, the relation of components entering 
into fluid mixtures made under conditions 
ranging from no efflux during the admission to 
that where efflux is equal to influx. For the latter 
condition, R=1, and Eq. (8) reduces to, 


y = —log. [(x2—x)/(x2—%1) J. (6) 


For no efflux, R=0, and Eq. (8) as well as its 
equivalent form reduces to the respective general 
equations given above for this condition of efflux. 

For any constant ratio of efflux to influx, if the 
rate at which the addition being made is known, y 
can be expressed as a function of the time ¢ 
elapsed after the addition began or, y=(f)t. 
However, irrespective-of what this relation may 
be, the concentration of a component at any 
time will be equal to that computed according to 
the methods outlined above as being due to a 
single instantaneous addition of the total 
quantity added up to the time, assuming as 
before that diffusion is complete and instan- 
taneous. 

If the ratio R of efflux to influx is not constant 
and its law of variation with respect to either y 
or ¢ can be stated, it may be possible to express 
the resulting concentration in terms of y or f, 
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depending on the form of the resulting ex- 
pressions. It will in general be found that with an 
equation thus derived, particularly if determi- 
nations of influx, efflux, or both are based on the 
theoretical laws of fluid flow, less accurate 
results are obtainable than by direct determi- 
nations thereof and application of Eqs. (3, 4, 8, 
or 8’) according to the method outlined in 
Section IT-3. 

Considering that the influx tends to impel a 
portion of the contents of the container ahead of 
it, the efflux is likely to contain less of the 
component being added than its average concen- 
tration at the given instant for the space as a 
whole. Hence for conditions of free or partial 
efflux, the average concentration of the added 
component may accordingly be greater than as 
computed according to the above development, 
which assumes complete and _ instantaneous 
diffusion. The relative position of inlets and 
outlets for mixture operations with given fluids 
can generally be arranged so that this condition 
will be assured. 


IV. DILUTION WHERE COMPONENT DILUTED Is 
Not PRESENT IN DILUENT 


While relations for dilution under mixture 
conditions ranging from free to no efflux can be 
obtained from the treatment in the preceding 
section, further consideration will be given to the 
most common case of dilution where the diluent 
contains none of the component that is to be 
diluted. For the case of free efflux Eq. (6) applies 
with x.=0. If for the sake of distinction C, be 
termed the initial concentration and C the 
concentration after dilution, Eq. (6) becomes, 


y= —log,. (C/C) = —2.303 logio (C/Ci). (9) 


This can be obtained readily also from the 
relation, 
dC=—CdyV/V, 
by integrating for y between the limits C and C,. 
Fig. 2 gives the relation between volumes of 
diluent and the reduced ratio C of the component 
after being diluted from initial concentration C). 
To extend the field of the diagram, two scales 
have been given for respective application to 
high and low initial concentrations of the 
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ConcenTaation C AFTER DILUTION LOW RANGE 


Fic. 2. Relation of volumes of diluting medium to concen- 
tration after dilution, assuming free efflux. 


component to be diluted. If necessary scales 
successively lower by a divisor of 10 can be 
assigned until the desired refinement is realized. 
Thus for the dilution of a component of concen- 
tration 0.60 to 0.0001, 2.30 volumes of diluent is 
indicated by the upper scale in Fig. 2 as needed to 
reduce its concentration to 0.06. With further 
addition of 2.30 volumes the concentration will 
be reduced to 0.006 (lower scale), and for a 
further addition of the same amount, to 0.0006, 
the latter being obtained with a one-tenth 
reduction of the scale. The reduction to 0.0001 
requires 1.8 additional volumes, which gives a 
total of 8.70 volumes of diluent for dilution from 
0.60 to 0.0001. 

For no efflux the relation can be obtained from 
the pertaining general equation given in the 
preceding section or directly from the obvious 
relation, 


C=C,/(1+y), whence 


vy=(Ci,—C)/C. 


For intermediate conditions of efflux Eq. (8) 
reduces to 


log.[1+y(1—R)]/(1—R) = —log,(C/C:) (10) 
and Eq. (8’) to' 
y=[1/(1—R)JE(C,/C)'-*-1]. = (10’) 


If R varies, successive computations, by using 
Eq. (10) or (10’), can be made for intervals over 
which R can be assumed as being constant, 
following the method outlined in Section II-3. 
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V. DILUTION OF OXYGEN OF THE AIR 


The above treatment can be conveniently 
applied where the oxygen in the normal air 
within a space is diluted such as by the admisston 
of inert or partially inert gas. By using the 
symbol z for oxygen concentration and by 
substituting in Eq. (6), 0.21 for x, as the initial 
concentration of oxygen in air, and x2=22, the 
concentration of oxygen in the diluting medium 
applied (flue gas or exhaust gas from internal 
combustion engines if used may contain some 
oxygen), the relation of quantity of diluent and 
oxygen concentration for the condition of free 
efflux is given by, 


v= —log, (——-) 
0.21 — Ze 


2— Ze 


u ). (11) 
0.21—25 


(11’) 








= — 2.303 logio ( 


This becomes 
y= —log, (2/2,) = 2.303 logio (2/21) 


if the dilution is made with a fully inert gas. 
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OXYGEN RATIO, Z 


Fic, 3. Relation of volumes of inert gas admitted to resulting 
oxygen ratio, assuming free efflux. 


Relations for Eq. (11) are given in Fig. 3 for 
values of ze of 0, 0.02, 0.05, and 0.10, concen- 
trations of oxygen in the diluting medium used. 

If no efflux takes place while gas is admitted, 
such as into a closed container, 
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5. Volume of air required to increase the oxygen ratio to value z from initial value 2; under 
conditions of free efflux. 
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Fic. 7. 

Fic. 6. Carbon dioxide concentration resulting from admitting one pound of carbon dioxide per 

; given volume units. ; 

Fic. 7, Oxygen ratio resulting from admitting one pound of carbon dioxide per given volume units, 
and corresponding carbon dioxide concentration. 


69 








70 S. H. 


z=(0.21+y22)/(y+1) or y=(0.21—2)/(¢—22), 


which becomes, y = (0.21 —2) /z, for dilution with 
a fully inert gas. 

These relations are given in Fig. 4 for assumed 
oxygen concentrations 22 of 0, 0.02, 0.05, and 0.10 
in the gas applied to reduce the normal oxygen 
concentration of air to the value z. 

For restoring the oxygen content after flooding 
with inert gas, the space may be ventilated by 
using normal air. Applying Eq. (6) with x2.=0.21, 
the number of air changes y for increase of 
oxygen content to z from initial values 2, will be, 


0.21—2 
= —_ : 
0.21 am Be 


0.21—2 
~ 2.303 logio| = - | (12) 
0.21—2, 


The curves in Fig. 5 give the ventilation 
requirements to obtain a desired oxygen ratio 2 
with initial oxygen ratios of 0, 0.05, 0.10, 0.15, 
and 0.20, respectively, assuming free efflux during 
the admission. 


VI. APPLICATIONS TO THE USE OF INERT GAS IN 
FIRE EXTINGUISHMENT 


Where a wholly inert gas such as carbon 
dioxide or nitrogen is used, the concentration 
under the condition of efflux present can be 
obtained from the developments in Section II 
and the curves in Fig. 1. While the condition of 
free efflux will generally be present, the corre- 
sponding values for an assumed condition of no 
efflux may be of interest. The resulting oxygen 
concentrations are obtained for conditions of free 
efflux from Eq. (11’) and Fig. 3. Comparative 
values for no efflux are given in Fig. 4. 

Where carbon dioxide stored in liquid form is 
used as the inert gas, the degree of protection 
obtained is conveniently expressed in terms of 
the number of cubic feet A of room volume 
allowed per pound of carbon dioxide provided 
as protection. If the room temperature after the 
drop due to the admission is assumed to rise 
again to a little over 30°C, the weight of the 
carbon dioxide gas at this temperature will be 
near 0.11 pound per cubic foot, or 9 cubic feet 
per pound. If y be the resulting number or 
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fractional number of volumes V of gas admitted 
each equal to the room volume, then 


yV/V=9/A. 


Substituting for y in Eq. (1), for the condition 
of free efflux, 


A = —9/[log. (1—x) ] 


The upper curve in Fig. 6 gives carbon dioxide 
concentrations corresponding to given room 
volumes allowed per pound of carbon dioxide 
applied as protection. 

For an assumed condition of no efflux ex- 
pressed by the equation, 


y=x (1—-x), 
the substitution y=9/A gives 
A=9[(1—x) /x], 


which relation is given by the lower curve in 
Fig. 6. 

The resulting reduced oxygen concentrations 
can be obtained by a similar substitution in Eq. 
(11’), whence for free efflux 


A=—9/[log. (2/0.21)] 
= 3.9086 /[logio (z, 0.21) ]. (14) 


This relation is given by the upper curve in 
Fig. 7 and that for no efflux by the lower curve, 
the equation for which is 


A =9[2/(0.21—2) ] 


obtained by substituting y =9/A in the applicable 
general equation, 


y = (0.21 —2)/z. 


For any condition of efflux, the relation be- 
tween oxygen concentration z and the added 
inert gas concentration x is expressed by, 


2=0.21(1—x), 


which is represented by the straight line in Fig. 7. 
After a space has been deluged with inert gas 
the ventilation requirements for increasing the 
oxygen content to a given limit making it safe 
for entry, can be obtained from Eq. (12) and 
Fig. 5. Similarly, by means of Eq. (5) and Fig. 2, 
the number of air changes can be computed that 
are needed to reduce known concentrations of 
inert or toxic gases to recognized safe limits. 
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The Flow of Compressible Fluids Through Porous Media and Some Problems in Heat 
Conduction 
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In Part I the general partial differential equations are 
derived for the flow of compressible liquids and of gases 
through porous media, on the basis of the generalized 
Darcy law that: for all homogeneous fluids the fluid 
velocity is directly proportional to the pressure gradient. 
The equation for compressible liquids turns out to be the 
Fourier heat conduction equation with the density as the 
dependent variable, while for gases a nonlinear parabolic 
equation is obtained. The steady state solutions for these 
equations for linear and plane radial flow are derived and 
discussed in Part II. The heat conduction equation, in 
polar coordinates with radial symmetry, for compressible 
liquids, is then solved in Part III for systems in which (a) 
the density is specified over both concentric circular 
boundaries of an annular region, (b) the density is given 
over one boundary and the normal derivative over the 
other, and (c) the normal derivative is given over both 
boundaries. For the last case a new elementary solution is 
introduced, not given in the standard English texts. The 
cases where the internal boundary reduces to an infini- 
tesimal sink are also given explicit solutions. Applications 


are made to (1) production history from a well whose 
pressure is reduced discontinuously to its final value; (2) 
production history from a well whose pressure is lowered so 
that the fluid density at the well drops linearly with the 
time; (3) the pressure rise in a well after shutting in; (4) 
pressure decline in the ‘East Texas’”’ oil field; (5) pro- 
duction decline from a single well at constant pressure in a 
closed reservoir; and (6) pressure decline in a closed reser- 
voir drained at a constant rate by a central well. The 
Green’s function, corresponding to a well displaced from 
the center of a closed reservoir, is derived and is applied to 
the problem of the interference between two wells draining 
the same circular reservoir. The solutions for all these 
problems correspond to similar heat conduction problems 
which have not been previously solved explicitly in the 
literature. In Part IV an approximation theory is given for 
the non-isothermal gas flow from a closed circular reservoir 
into a central well, and is illustrated by a numerical example 
showing the pressure and production decline for such a 
reservoir on a sudden lowering of the well pressure from an 
initially uniform value. 





INTRODUCTION: FUNDAMENTAL EQUATIONS 


I. 
IN several previous papers from the Gulf 
Research and Development Corporation in 
this Journal and in the A.J.M.M.E. Transactions 
treatments, experimental and theoretical, have 
been given of various problems involving the 
flow through sands or porous media of ‘‘dead”’ 
incompressible liquids or of ideal gases, under 
isothermal conditions. The solutions for these 
problems have been developed not only because 
they represented approximations to situations 
which actually arise in practical petroleum 
technology, but also because they form an ideal- 
ized framework for at least the qualitative 
discussion of the vastly more difficult problems 
involving the flow of gas oil mixtures, thé 
understanding of which is really the essential 
purpose of the study of the hydrodynamics of the 
flow of fluids through porous media as related to 
the oil industry. 


For the same reasons it is of interest to examine 
in detail the properties of flow systems in which 
the porous media carry general compressible 
fluids, such as ideal gases in non-isothermal 
expansion, and liquids having nonvanishing 
compressibilities. With regard to their immediate 
applicability, the case of non-isothermal ex- 
pansion of ideal gases is involved in the pro- 
duction, from gas wells which flow at high rates, 
while the liquid compressibility must be con- 
sidered in those problems where the oil field is 
being produced by a water drive, for it is clear 
that the compressibility of the water in the sands 
adjacent to the oil sands will provide the source 
for at least part of the energy driving the oil into 
the producing wells. But even aside from these 
specific applications the analysis of the flow of 
these fluids through porous media is of value in 
that they give the nature of the flow for the 
extreme and idealized types of fluids which are 








evidently the component parts of the gas-oil 
mixture fluids occurring in practice. And in fact 
the ideal compressible liquid which already 
contains characteristics of both the gas and liquid 
systems should really give at least a qualitative 
representation of some of the new features to be 
expected of the gas-liquid mixture fluid. 

In developing a theory of the flow of com- 
pressible fluids through porous media it is 
convenient to follow the well-known methods and 
principles of the classical hydrodynamics. In 
fact, the flow of fluids through porous media 
should, in principle, be describable by the laws of 
classical hydrodynamics and hence require no 
new formulation or treatment. It is, however, 
only because of the hopelessly insurmountable 
analytical difficulties of the solution of the 
Stokes-Navier equations for such complex geo- 
metrical systems as are formed by the tortuous 
paths in a body of sand that a macroscopic theory 
which averages out these individual fluid chan- 
nels, and is based on direct empirical laws derived 
from experiments with sand bodies, must be used 
in order to make any progress in the study of 
porous media as carriers of viscous fluids. 

This macroscopic theory, on the other hand, 
differs from the Stokes-Navier equations only in 
the dynamical characterization of the hydro- 
dynamics. The equation of continuity and the 
equation of state defining the fluid must remain 
valid, of course, regardless of the details of the 
hydrodynamics. Thus, if y, 6 are the fluid 
density and vector velocity at (x, y, z) and at the 
time ¢, the equation of continuity requires that: 


div. (yi) = —fdy/dt, (1) 


where f is the porosity—fractional open space—of 
the medium. Further, the equation of state 
defining the particular fluid of interest together 
with the equation describing the thermodynamic 
character of the flow will give a relation of the 
form: 


y=y(b), (2) 


where ? is the pressure at (x, y, 2) and at ¢. 

With regard to the dynamical equivalent of the 
Stokes-Navier equations our fundamental law 
will be given by: 
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b= —(k/u)Vp, (3)! 


where u is the viscosity of the fluid and is the 
“‘permeability’’ of the porous medium, and is 
perhaps best defined by Eq. (3) itself. In the 
case of the liquids, this law was established as 
long ago as 1856? by H. D’Arcy and shown to 
hold quite accurately as long as the velocities of 
flow did not become excessive and turbulence 
did not set in. It has been repeatedly confirmed 
since then by other investigators. For gases, Eq. 
(3) has again been empirically confirmed, al- 
though definite proof as to its validity has been 
obtained only recently.* Here too, deviations 
were found, even in the case of consolidated 
sands, but only at such high velocities as would 
be found to occur only in wells producing gas at 
very high rates, and then the deviations would be 
confined to the immediate vicinities of the well 
bores. , 

Although no experiments seem to have ever 
been made with liquids to test Eq. (3) over 
pressure ranges in which the compressibility 
effects would become appreciable, it seems safe to 
assume, that since it does hold for the extreme 
cases of gases and effectively incompressible 
liquids, it also holds for liquids of high com- 
pressibilities, or over pressure ranges in which the 
compressibility would have appreciable effects 
upon the fluid density. It shall, therefore, be 
supposed that Eq. (3) holds for all homogeneous 
fluids, flowing under viscous conditions—absence 
of kinetic energy and turbulence effects. 

There is now left simply the specification of the 
nature of the fluid and the thermodynamic 
character of the flow—Egq. (2)—to obtain a 
differential equation in the single variables y or p. 
For with Eq. (3), Eq. (1) becomes: 


div. (ykVp) =fudy/at. (4) 


As a sufficiently general equation including all 
fluids and conditions of flow of practical interest 
one may take: 


. | When gravity is taken into account p must be replaced 
by the function p—~ygz, where g is the acceleration of 
gravity and z the vertical coordinate, the — sign implying 
the downward direction of +z. ’ 

2H. D’Arcy, Les Fontaines Publiques dela Ville de Dijon, 
1856. 

*M. Muskat and H. G. Botset, Physics 1, 27 (1931). 
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Y= yopme’”. (5) 


The particular fluids of physical significance may 
now be classified as follows: 


liquids: m=0 
incompressible liquids: 8 =0 
compressible liquids: 8 +0* 

gases: B=0 
isothermal expansion: m= 1 


’ : : specific heat at const. vol. 
adiabatic expansion: m= 





specific heat at const. pres.. 


It follows at once that for strictly incompressible 
liquids there can be no time transients, or non- 
steady states, within the system unless the 
boundary conditions vary with the time. For 
homogeneous sands (k=const.) the pressure 
obeys Laplace’s equation. For all other cases the 
non-steady state systems are governed by second 
order parabolic partial differential equations, 
which reduce to Laplace’s equation for the steady 
state. Furthermore, for compressible liquids this 
equation, when using the density as the depend- 
ent variable, is identical with that for heat 
conduction. Thus, for compressible liquids one 
obtains: 


V7 = (fBu/k)(dy dt). (6) 


For gases on the other hand, Eq. (4) takes the 
form of a nonlinear equation as: 


(1+m)pyo"!”™ f dy 


Vy tm) /m = ——___a—-—____ oe (7) 


’ 


k al 

both cases referring to homogeneous sands. 

Before entering upon a detailed analysis of 
these equations it may be of interest to compare 
in a general way the rapidity of the transients in 
gas and compressible liquid systems. By taking 
the particular case of isothermal gas flow (m=1), 
it is seen from Eq. (7) that the time will enter in 
the expression for 7? essentially as kt/fyou. Now 
for a given sand f and k will be independent of the 
fluid and for ideal gases yp = w/RT where R is the 
gas constant per mol, w the molecular weight and 


‘ Although true liquids show a compressibility decreasing 
with the pressure its variation in the pressure range of 
practical interest (1-200 atmos.) is rather small and will be 
neglected. 


T the absolute temperature. Hence, for isother- 
mal gas flow ¢ will enter essentially as tRT wy, so 
that the time rate of change in the system will be 
proportional to RT/wy. For compressible liquids, 
on the other hand, it follows similarly from Eq. 
(6) that the time rate of change will be pro- 
portional to 1/p 8. 

Now for a gas such as methane R7T/w~ 1500, 
whereas 8 has a value for mineral oils of the order 
of 10-4. The viscosity, however, for mineral oils is 
of the order of 10° times that of methane. Hence 
from the resultant of these two factors it follows 
that the rate of decline in a gas reservoir under 
isothermal conditions will be of the order of 100 
times as large as in one containing a gas-free but 
compressible oil. Apparently then, the low 
viscosity of gases more than counterbalances 
their high elasticity in giving rise to very rapid 
transients.® 


II. StEADY STATE FLOW 


When the flow systems are maintained in a 
steady state, the right sides of Eqs. (6) and (7) 
may be set equal to zero, thus giving: 


v-y7=0 (6a) 
and 
Vyiitm m= (), (7a) 


Hence y, in the case of compressible liquids, 
and y“+™)/™ in the case of gases, take the réle of 
potential functions for steady state systems and 
may be treated by the usual methods of potential 
theory. When 7 is found the pressure distribution 
may be computed by returning to Eq. (5) and 
solving for p. 

Thus in the case of linear systems it is readily 
found that if 71, y2 and pi, pe are the densities 
and corresponding pressures at the extremeties 
of a linear flow channel of length LZ, the distri- 
butions in y and p are given by: 


°> The fact that Eq. (7) involves 7? (for m=1) whereas 
Eq. (6) is linear in y does not invalidate this comparison. 
For if Eq. (5) is approximated by y=y0(1+8p) it is found 
that Eq. (4) takes exactly the same form for compressible 
liquids-as for gases, i.e., it also involves V?7? rather than 
Vy. 
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Compressible liquids: 


y= —Llv2—-71)/L ]x+72, (8) 
p=(1 8) log |}—[(yve-y1) L Jxt+ye2} —(1 8B) log yo. (8a) 
Gases: 
¥ 1+m)/m — —((y20*™ m — a (l+m) m)/L | x+vye l+m m (9) 
gota = —[(p'*™— p,'**) L | x+p.'*™, (9a) 


where x is measured from the point where y and p have the values 72, ps. 
The rates of mass flow for these systems are given by the equation: 


(7- (k (wB)(dy Ax): compressible liquids 
Q= —(ky, u)(Ap/dx) (10) 
= —(k/(1+m)pyo!'”)(0yr™!™ Ax): gases! 
per unit cross section, k being the permeability of the medium and x the viscosity of the fluid. 
By applying this to Eqs. (8)—(9), it is found that for: 
Compressible liquids: 
Q=(k, Bu) (vy2—¥1)/ L= (Ryo BuL)(e8?2 — e891) = Qo(1+ Bp) (11) 


to terms in 8°, where Q is the flow for the case of an incompressible liquid of density yo, and viscosity 
u, and pis the algebraic mean pressure. 


Gases: 


R(y_"tt™ m — oy (lm) m kyo ? ; 1—(p, p2)'*™ 
Q= wetie - bit se SOOT etna aimee (12) 
(1+m)pyol™L w(l+m)L 1+m 1—(p; ps)? 








where Qp is the isothermal flow rate (m =1), after observing that y)~p.~” if the densities for the two 
values of m are to be the same at the pressure po. 

Similarly, in the case of radial flow systems, if y.., 6. denote the density and pressure at the well 
of radius r., and y., p. the corresponding values at a distant radius (from the center of the well) r., 
the density and pressure at any distance 7, are given by: 


Compressible liquids: 








e 
¥ =D —Yw) log 7/ re | log re Tut Ves (13) 
1 f vole?” —erw) log rr, 1 
p=-— log ;— +747 —— log yo. (13a) 
B | log ry 8 
Gases: 
(>, “it™ m— my, (1+m) m) log |e 
it~ m — +y¥,(°™ me (14) 
log rr, 
(p3t™ — p,'*™) log v/v. 
pit" =—_—__—___ — ——+p,!*™, (14a) 





log fr. x 
The mass outflow rate per unit sand thickness corresponding to these distributions is given by: 


2rrk dv 
= —: compressible liquids 
2rrky ap uB or | 

Q = at > (15) 
um or 2rrk oy !rmi= 
=——______ ———-: gases | 
(1-++-m)pyo''™ or / 
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OF 
which leads to: 


Compressible liquids: 


Q= (2k / UB) (Ye — Yw) / (Og Fe/ Pw) =Qo(1 + Bp) 


COMPRESSIBLE 





FLUIDS 


~I 
cn 


(16) 


to terms in 6%, Qy being the radial flux for an incompressible liquid. 


Gases: 


Q=(2Q0/(1+m))[1 — (pu/pe)'*"]/(1 — (pu / pe)? J, 


(17) 


where again Q, is the isothermal gas flow rate (m=1). 


Since yo(1+  ) is, to terms in 6”, the average 
density of the liquid in the system, Eqs. (11) and 
(16) show that the flux in a steady state com- 
pressible liquid system may be computed, to 
terms in 6°, by taking for the density the 
algebraic mean of the densities at the two 
extremeties of the system but otherwise using 
the ordinary formulae where the compressibility 
is neglected. 

In the case of gases, however, there is a definite 
effect of the deviation from isothermal flow 
(m=1). For as is readily seen from Eqs. (12) 
and (17), Q, for non-isothermal flow (m<1), is 
greater than Qo, the isothermal flow. Moreover, 
Q/Qo increases as the deviation increases, i.e., as 
m decreases from the value 1. Furthermore, for a 
given type of flow, fixed m, Q/Qo increases as 
Pw/ pe decreases—as the pressure differential and 
production rates increase. 

The actual order of magnitude of the effect of 
the compressibility is quite small. Thus, taking a 
compressibility of 5X 10~‘ atmos. (attributing a 
factor of 10 to the effect of dispersed gases) and 
an average pressure of 100 atmos. the factor 
1+) becomes 1.05, or an effect of only 5 percent. 
Similarly the pressure distributions correspond- 
ing to Eqs. (8a) and (13a) are almost indis- 
tinguishable from those for incompressible liquid 
systems, except for large 8p, even though the 
equations for the latter systems are much 
simpler than Eqs. (8a) and (13a). 

As an example for the case of gases one may 
take m=0.71 which, for air, corresponds to 
adiabatic expansion. Then if p,,/p.=0.1 it is 
found that Q/Qo=1.16. Thus, there will be 16 
percent more fluid passing through the system if 
it expands adiabatically (and p,./p.=0.1) than 
if the expansion is isothermal. It is to be noted 





that the main reason for this is the fact that the 
density at the outlet is higher, for equal inflow 
densities, for the non-isothermal flow than for the 
case of isothermal flow. This excess more than 
compensates for the lower outlet pressure 
gradients of the non-isothermal flow. This may be 
seen from the pressure distribution curves plotted 
in Fig. 1 for a linear system for the adiabatic case 
of m=0.71 and for the isothermal case of m=1. 














Fic. 1. Pressure distribution in a linear channel flowing 
gas under isothermal conditions (curve I) and adiabatic 
conditions (for air, curve II). x/L=fractional distance 
along channel; p/p2= (pressure at x/L)/(pressure at x =0). 


III. Non-STEADY STATES: COMPRESSIBLE LIQUIDS 


When the flow systems are in transient con- 
ditions the general Eqs. (6) and (7) must be used. 
In the case of compressible liquids, this requires 
the treatment of the heat conduction Eq. (6) 
which can be accomplished by an application and 
extension of some well-known theory. Eq. (7), on 
the other hand, is nonlinear and does not permit 
an explicit solution. An approximation method 
will, therefore, be developed, which is analogous 
to that already given for the case of isothermal 
flow, in a previous paper.* Because of this 
fundamental difference, however, it will be 
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convenient to treat these problems separately, 
and as the heat conduction equation can be 
solved rigorously it will be considered first. In 
each case, unless otherwise noted, the treatment 
will be confined to two dimensional systems with 
radial symmetry, as these are of most practical 
interest. 


A. Types of problems; fundamental solutions 


There are three types of problems of practical 
interest which may be proposed with respect to 
radial flow systems. They may be stated ana- 
lytically in the forms: 


(1): y=filt) :r=nry 
= fo(t) :r=fre (18) 

=gir) :i=0, 

(2): y=filt) :r=n 
r dy, Or=f2(t) :r=Pe (19) 

y=x2(r) :1=0, 

(3): ray dr=fi(t) :r=nri 
= felt) :r=Pre (20) 

y=gz(r) : t=0, 
where 7r;, 72 define two concentric circular 


boundaries limiting the region of interest, ¢ is the 
time variable, and f), fe, g are arbitrary functions 
to be preassigned so as to correspond to the 
particular physical problem under consideration. 

Since the densities y uniquely determine the 
pressure p, the above equations may be directly 
interpreted in terms of pressure distributions. 
Thus, case (1) corresponds to the situation in 
which the pressures at the two boundaries 7 =1r;, 
r=re—one of which may be chosen as a well 
radius—are made to vary (or kept constant) in a 
preassigned manner, the sand having a known, 
but arbitrary pressure distribution at the initial 
instant. The solution to this problem would yield 
not only the pressure distribution at any later 
time but also the fluxes entering or leaving the 
boundaries of the system. 

Case (2) involves a situation in which the flux 
rate at one of the boundaries is known and the 
pressure at the other. As an example one may 
suppose that the flux rate at the external 
boundary is zero so that the sand forms a closed 
system. The internal radius might be taken to 
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define a well at which the pressure is varied (or 
held constant) in a_ specified manner. The 
solution of the problem would then serve to give 
the production rate from the well under these 
conditions. On the other hand, the pressure 
might be preassigned at the external boundary 
and the flux rate at the internal boundary, the 
solution showing how the pressure will vary at the 
internal boundary. Such an interpretation may 
be applied—as will be seen later—to the pressure 
decline in a field as the ‘‘East Texas”’ field on the 
basis of certain assumptions as to the basic 
mechanism of the production. 

Case (3), finally, corresponds to a system in 
which the flux rate is specified at both boundaries. 
For example, the flux at the external boundary 
may be taken as zero—a closed system—and 
the production rate from the well which may be 
assumed as the inner radius, may be arbitrarily 
specified. The pressure decline at the well, will 
then be given by the solution to the problem of 
Case (3). 

These solutions must, of course, all be derived 
from the fundamental differential Eq. (6) which 
may be written, in polar coordinates as: 
Oy 1 ay) Oy k 

-— e——: gm 


=——, (21) 
rar) at ufB 


il ein 
or 

Although the heat conduction equation has 
been given an extensive treatment in the 
literature, the solutions for the particular 
problems of interest here do not seem to have 
been developed in detail. Thus, Carslaw® gives 
the explicit solution for Eq. (18) only if f;=fe=0 
and the special case of Eq. (19) where 7,;=0, 
f2=0, and makes no direct mention at all of 
solutions for the problems of Eq. (20). However, 
since the solutions to be presented here may be 
obtained by extending the methods developed 
by Carslaw only the essential intermediate steps 
of the derivations will be given below. 

Before proceeding now with the detailed 
solutions for the above problems it will be 
convenient to summarize some analytical pre- 
liminaries which will be extensively used in the 
following developments. First it is to be observed 


°H. S. Carslaw, Mathematical Theory of the Conduction 
of Heat in Solids, 2nd Ed. 1921. 
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that we have available for the treatment of the 


problems of Eqs. (18)—(20) the following three 
elementary solutions: 


7 = log r, (22) 

y=r+4ut, (23)? 
Jv\(ar) 

y= oe, (24) 
Y,(ar) 


where @ is any real number and Jo, Yo are the . 


zero order Bessel functions of the first and second 
kind. 

Combining these elementary solutions to 
correspond to problems with boundary values 
independent of the.time the result will be used 
that if yo(r, 4) is a solution of Eq. (21) which 
vanishes at ¢=0 then the function: 


t 


0 





a ave—a 


b 1; dU i 
f rer(arydr= a PU*ar) + —(an)| . 
ar a 


a a 


The explicit solutions will now be presented: 


B. Solution for Case (1) 
Let solutions y1, y2 and y; be defined so that: 


(v1, Y2, Ys) = (1, 0, 0) 
= (0, 1, 0) 
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is also a solution of Eq. (21). This is a special 

case of Duhamel’s theorem and may be readily 

verified by substitution. This formula permits the 

derivation of a solution which corresponds to 

variable boundary conditions from one in which 

the boundary values are constants. 
Now let: 

U(ar) =AJ (ar) +BY,(ar), (26) 
where A and B are any arbitrary constants which 
may depend on a, but not on r. Then the 
following integrals may be easily verified: 


b 


f rU(ar)dr=(—1 a)[rdU/dr],°, (27) 


a 


ab 1 dU P 
| rU(ar) log rd =| U-r tog | , (28) 
” a dr " 


9 
2 


ad 1 dU dU > 
f rU(ar) U(a'r)dr= —— | Utara) _ rU(a'n)—(ar) | . axa’, 


dr dr 


(29) 


(30) 


~-r=?r; 


:7=T7o>. 


(31) 


= (0, 0, g(r)) : t=0 


Then Eq. (18) will be satisfied by: 


; Oy1 Ove 
Y=73(", n+ f | fun), t—d)+ fe(A)—(r, t—X) dX. 
0 ot ot 


Now let: 


U(a,r) = Yo (@n%2)J olan?) mad J (p12) Y,(a,r), 


where a,, is a root of: 


U(a,r1) =0. 


7 It may be of interest to note that the simple solution of 
Eq. (23) is not even mentioned in the standard (in English) 
texts on heat conduction, although it is necessary for the 
solution of the problems of Eq. (20). Similarly neither the 
corresponding solution: y =(1/r)(7°+2xt) for systems with 


(32) 


(33) 


(33a) 





spherical symmetry which would be required for the 
spherical problem analogous to Eq. (20), nor the Cartesian 
coordinate solutions 5,"q.?+2nxt have been found by the 
author in the standard texts. 
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Applying Eqs. (27)—(30), noting that: 
U"(a,%2) = —2/are : U' (ants) = — 2S o(anr2)/ HJ o(ants) (33b) 
and expanding 7; formally as: 


v1 = (log 7/72) / (log r1/r2) +A, Ula, rye *e""* (34)8 
it is found that: ; 
log r/ re Slant) J o(anre) Ula,rye—*e"! 
ae ke > 7 - (34a) 
log Y,/ Te JP (ant1) — J 7 (an%2) 





Expanding v2 similarly it takes the form: 


log 7/7) J (ant) U(a,r)e**""* 
fg Wh ce emareemnsnen eae . (35) 
log 7; /f2 J? (anti) — J?(a,12) 





v3, finally, may be expanded in a.similar manner with the result: 


Tan? J (ants) Ul(agrye*"* >? 
¥s=—L : : { re(r)Ula,r)dr. (36) 
2 J Pants) — J? (anP2) 1 





Adding these solutions together, as is indicated in Eq. (32), it is found that: 
a ned ol Q@n%1) Ul(a,rje seit ios 


Y="L- - 
JP (ant) — Je? (anPe) L2 





3 
Jalen) { rg(r) U(a,r)dr 
"1 


t t (37) 
= 1Jolacrs) furyew?ar aJolasts) [ fainvereetn} 
0 0 


It should be observed, in using this equation, that to check for the boundary values one must not 
replace r directly by r=r1, 72, for then in virtue of Eqs. (33), (33a) y will vanish identically. One 
should rather sum the series first and then make r approach the limiting values of 7; or 72. This 
apparent difficulty arises from the fact that in setting up the expressions of Eqs. (34a) and (35) for 
y; and y2 expansions of the form: 

log r/ri, 2= 0 A,U (aur) 


were used, for which the series at the right are discontinuous at r=/7e, 7;. However, if the series be 
summed first and the boundary values substituted in the sums no difficulties will arise. 


(1) Example of Case (1) 

As an illustration of Eq. (37) one may consider the problem of an oil field produced by a water 
drive, the water extending over a large area outside the oil pool. For simplicity the system of wells 
drilled into the oil field will be replaced by a single circular well of radius equal to a dimension of the 
field. The initial density distribution g(r), will be taken as constant. The radius 72 will be considered as 
the external radius r,, the pressure at which—and hence the density—will be taken as constant, and 
equal to the initial pressure. It will be furthermore supposed that at the boundary of the field, defined 
by the ‘‘well radius”’ 7,,(7,), the pressure in one case, has been suddenly dropped to ~,..(y =7.), and in 
the other case that the pressure declines continuously so that the density falls linearly to yo. (Cf. 
Fig. 2.) 

Analytically these assumptions imply that: 


8 It should be observed that in the series of Eq. (34) and in all the corresponding series to follow, the summations 
are supposed to extend only over the positive roots a». 
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y=g(r)=yis (p=pi) :t=0, ) 
Y= ¥e = Vii (pe= pi) :r=rA(r2), 
; } (38) 
(a): Y=%es (Pw = (1/8) log +. Yo) 
r=7,(f)) 
(b): yY=7i--el; (p,.=(1 B) log [(yi-—el) vo ]) “4 


By putting these into Eq. (37) and summing the series not containing the time it is found that: 


vi log 7/TwtYu log r./r J (cente) I Oe ntw) U (annem *2n*! 
a): y= annie) \i ia ., ews —., (39) 
log Le, Tw J Pant) aa J F(a nPe) 








yi log r/ftut(yi-—et) logr./r me Solent.) Jo(anfw) U(a,7) 


( b) ° ee ee 








log 7../ Tw Kae | IJ (anP i) — J 0? (Ante) } 
te J y(@n%)Jo(@n1,-) U (annem "2"! 
+—2, - (40) 
kK an? | JPR (anfw) —J 0 (Qnte) | 


That these equations satisfy the boundary and initial conditions of Eq. (38) may be verified by 
direct substitution. 


Of particular practical interest is the flux or production rate that would be obtained under the above 
conditions. This will evidently be given by: 
Q=(2Qrk/p)(ry 0P/87r) re =2TR/wB re (OY /OP) rw. (41) 


By applying this to Eqs. (39) and (40), summing, and introducing the dimensionless notation: 


Te/ T=; Ono =Xn3 xt/ty2=t (42) 
it is found that: 








2ak(yi— Vw) 1 J2(x,p)e-2""* 
-+2>° , (43) 


(a): Q=- x 
up Llog p J? (¥n) —JeP(xnp) 








t p” 1 1 1 log p\? JeP(xnp)e-2"* 
(b): Q= defer] —— + — —-————— log p— ( ) —2> | (44) 
log p 2(log p)*\2 2p? p? p Xn?) TP (xn) —J0?(xXnp) } 


In evaluating Eqs. (43) and (44) one must determine first the roots a, of Eq. (33a). In the notation 
of Eq. (42) the first seven roots for the case p=5(7,,~ 20 miles; r-~100 miles) are: 


n: 1 2 3 4 5 6 7 
Za: 0.7632 1.5575 2.3470 3.1352 3.9210 4.7073 5.4933 


The values of the Bessel functions may then be read off from the tables of Hayashi. The final results 
are given in Fig. 3 where the quantities QBu/2rk(y:i—yw) and Q/2zfer,2, denoted by Q, are plotted 
against ¢ as curves I and II, respectively. 

As should be expected, curve I begins with infinitely large values of Q and asymptotically ap- 
proaches the constant 1/log 5=0.62. This is the steady state flux that will pass through the system 
after the transient effects represented by the series in Eq. (43), due to the compressibility, have 
disappeared. Taking k/u=2 and B so small that yi—yu~708(Pe— Pw), the value Q=1, for curve I 
corresponds to about 200 bbl./day/ft. of sand per atmos. pressure differential (p.—p,=1 atmos). 
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t 
Fic. 2. Diagrammatic representation of an oil field pro- Fic. 3. Production history of fields produced by a water 
duced by a radial water drive. drive maintained at constant pressure. Curve I: Field 


pressure suddenly dropped to and maintained at p.; Curve 
II: Field pressure lowered continuously so that the fluid 
density falls linearly with the time. ¢ is the “dimensionless 
time” and Q the “dimensionless flux rates.” 


To get the numerical equivalent of the ¢ scale, values of the compressibility 8 and porosity f must 
be chosen. Taking f=0.2 and 8~4.5 X10-* atmos. (10 times that of water) one has: t~ 1.85 x 10-4 
(days) so that :=0.1 corresponds to approximately 1} years. This system will, therefore, possess 
very long period transients, for even after 6 years the flux rate would still be twice that in the steady 
state. 

In the case of curve II one unit on the scale of ordinates has the numerical equivalent of 4.06¢ 
x 104, ft. sand, or approximately 153,000 bbl. day /ft. sand for a pressure decline at 7, of 2 lb. day. 
However, if the initial pressure is 1600 lb. the rate of decline could be maintained for only 800 days, 
or to t~0.15 at which Q~0.5. Hence the production rate by the time the field pressure had been 
brought to zero would have risen to 75,000 bbl. /day ft. of sand. At the same time 0 ~ 1.9 for curve I, 
so that if 4p~ 1600 Ib. the system in which the pressure had been suddenly dropped to zero would be 
producing at the rate of 41,000 bbl. ‘day. 

With regard to the cumulative production it is of interest to note that: 


Ay Ap Ap 
Pia)= { Qla)dt=—Q0s) 2— —Q b)= —) h). 
0 € dp dt) 2 


de 
at 


where Q(a) and Q(d) represent the expressions of Eqs. (43) and (44), respectively. Hence the cumu- 
lative production for case (a) at #=0.15 is 


P(a) = 800 Q(5) =6 X10 bbl. ft. sand. 


P(b), on the other hand, may be most readily obtained from the area under curve II, the result 
being 4 x 10’ bbl. /ft. of sand. The general variation of P(b) with time will be essentially of the form of 
f, since after an initial short period transient Q(>) increases approximately linearly with ¢. For the 
same reason, P(a), in virtue of Eq. (45), will increase approximately linearly in ¢ after the initial 
transients have passed. 


(2) Limiting case of vanishing internal radius 


Although Eq. (37) is independent of the absolute values of 7,. and r,, one may not proceed to the 
limit p~ x, 7,.~0. For if the flux should be nonvanishing with a well of vanishing radius the well must 
he replaced by a mathematical sink, with a negative infinite pressure. The specification of finite fixed 











OW tte TF eee 


le 


It 
of 
1e 
al 


ne 


od 








FLOW OF COMPRESSIBLE FLUIDS 81 


or variable pressures becomes meaningless as r,,—0. If, however, one does not require a nonvanishing 
flux, the pressure or density at the well (7,,=0) becomes uniquely determined by that on the external 
boundary. The distribution in y is then given by: 


2 J (arr) ‘ *r, at , 
y=—--— e ~t | re(r)Jyla,r)dr+ wren (ast) f fainyernran (46) 
( 0 


r= JP(anr) ) 


where Jo(a,r-) =0, fo(t) is the value of y maintained at the external boundary r=r,, and g(r) is the 
initial distribution. 

An application of Eq. (46) of some practical interest may be found in the problem of taking 
“bottom hole’’ pressures in a producing oil field. The procedure is to send down to the bottom of a 
producing well a maximum reading bottom hole pressure gauge, close in the well and after a certain 
period remove the gauge and record the reading. This value is taken as the ‘‘reservoir’’ pressure in the 
field about the well chosen, and the question arises as to the time that is required after closing in for 
the pressure at the well to rise to a value appropriate to its “vicinity.” 

The following solution to this problem will not strictly apply to cases in which there is free gas in 
the sand about the well, due to the variation in the effective permeability of the sand as the pressure 
builds up and the gas is compressed. However, it will be valid in fields producing entirely by water 
drive—as seems to be the case in “‘East Texas’’—in which no gas is evolved from the oil until it 
reaches the well bore; and it will give a qualitative description of the pressure rise even in sands 
flowing both gas and oil, the curve of Fig. 4 falling more steeply in these cases, due to the rise in the 
effective permeability as the pressure about the well rises. 

Supposing then that the well, before closing in, has been flowing in a steady state condition, its 
density distribution at the initial instant, g(r) in Eq. (46) will be given by: 


=[(¥e—Ywi) log r/r./log 1. / tue lAYwis (47) 
where y, is the “reservoir density”’ to be determined, and y,,; is the flowing bottom hole density just 


before the test is made, 7, being the well radius. And since the field as a whole continues flowing f2(t) 
may be given the value 7. 


Under these conditions Eq. (46) takes the form: 
—[2(¥e— Yi) /12 log re, tw JECT ol(anren "08 / a2 Fant) |. (48) 


The pressure differential between the well (at r=r,,) and the radius r, is then given by: 


Pi—Pe Yew mt 
=F (49) 
Pp. — Pwi Ye Yui log p Xn” 2J “(z.) 





in the notation of Eq. (42), with p,, p. corresponding to y,, Yw. Thus the fractional rise in the pressure 
differential at any time is independent of the absolute values of the pressures and is determined only 
by the dimensionless constants in the series of Eq. (49). Eq. (49) is plotted in Fig. 4 for p= 2640, and 
rp for the factor 1/log p is a universal curve. Thus for k/u=0.2, f=0.2, 8=4.5 X10 atmos.~', 

= 660’, ¢=0.1~t=0.5 hours, so that 90 percent of the original pressure differential would be 
pita ed in about 3 hour and 99 percent in ‘23 hours. These values are of the same order of magnitude 


as those observed in the ‘‘East Texas’ field where the constants correspond approximately to those 
chosen above. 





— 
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C. Solution for Case (2)° 


For the problem of Eq. (19) it is again convenient to split up the final solution into three simple 
ones, 71, Y2 and y;3 where: 


(v1, ¥2, ¥3) = (1, 0, O) : r=r7; 
ro(d71/ Or, Ov2 Or, Oys/ dr) =(0, 1,0): ent (50) 
(v1, Y2, Ya) = (0, 0, g(r)) : t=0 J 
Then, similarly to Eq. (32), the final solution satisfying Eq. (19) may be written as: 
‘ Oy1(r, t—X) Oy2(r, t—dX) 
y¥=y¥3(r, t) +f E d) eet 7... see =|a. (51) 
0 ot or 
Introducing now the new function: 
U(a,r) = Yila,r2)Jo(anr) — Ji (are) Yo(a,7) (52) 


and choosing a, such that 
U(a,r:) =0 (32b) 


and then applying again Eqs. (27) to (30) the values of y1, v2, ys may at once be written, down as: 


Tolar Ii Q,72) U(a,r)e Kant 





n=i+e2— SP eereneeer (53) 
J ?(an0%1) —J7(a,%2) 
re. J 7(an%1) U(aprje*e""' 
¥2=log —+—>-— Cama a , (34) 
ri fe Gn} Je (ant)—JSP(anre) | 
rr? an? J (a7) U(anrye Kant ry 
¥3=—>)—_—_—-_—_ —_ —f rg(r)U(a,r)dr. (33) 
2 J (a,%1) —J7(a,%e) r 


i 


The resultant solution is then, by Eq. (51): 


t 


an? ST o(ant:) Ula,rye Kan®t us ls 9 
y=TL pares Jolauns) f re(r) U(anr)dr— KJ \(a,%2) f furyeres*nan 
J? (anti) —J7(anP2) 2 Tr ev 


- 1 





— kJ o(a@nP1); (aur) { fireman (56) 


0 


As illustrations of Eq. (56) treatments will be given first of a problem equivalent to that of the 
pressure decline in the ‘‘East Texas’’ field and then the pressure decline in a single well. 


(1) Pressure decline in the ‘East Texas”’ field 


It will be supposed that the “East Texas’’ field is being produced by a water drive and the geo- 
metrical conditions will be assumed to be the same as used in Case A, the system being represented 


® The special case in which the reservoir pressure and the 
well production rate are fixed with the solution given by for the results quoted in the work of Moore, Schilthuis 
Eq. (54) corresponds to the “‘constant production rate’’ and Hurst. Although the method of treatment is the 
case of Moore, Schilthuis and Hurst (Oil Weekly, May 22, same as that used here, the two problems solved by Hurst 
1933). Note added in proof: Since this paper was submitted, are both special cases of Eq. (56), one being given, except 
Hurst has published—Physics 5, 20 (1934)—the analysis for constants, by Eq. (53) and the other by Eq. (54). 
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by Fig. 2. It will be further supposed that the external pressure is kept fixed and equal to its initial 
value. The flux at the “‘well’’ or the production rate of the field will be taken as a} up to t=), b to 
t=te, c to t=ts, d to t=t, and e for ¢>¢,. Analytically, this implies: 
y=2(r)=7:; (p=pi) :t= 
Y=7e=7i3 (p.=pi) : r=r.(r1), 
a :0St=t, : 3/20/3134/1/33 
b 2 tStSte : 4/1/3354/22/33 
OY 
re) =4¢ :teStSts : 4/22/33-7/1/33+-r=r,(re). (; 
"|: tgStSty 2 7/1/33-39/1/33 


wi 
~I 
— 








le :tsSt :9/1/334 


4 


The first two conditions applied to Eq. (56) immediately give for the density difference between the 


two boundaries: 

2k JP? (anr.e~*™ 

n-ne =—E- —f fa(d)er@?rd (58) 
Tuo* Je? ( ante) — Ji? (Qn? w 





By using again the notation of Eq. (44) and introducing the function: 


7 J2(xnp)em=n"= | 
y(t) =2>° ' (59) 
ay } JP (Xnp) a J?(xn) ; 





the last of the conditions of Eq. (57) gives: 


0=t=h, Ay =a log p+ay(t); 
(60) 


hst=h, Ay =b log p+ay(t)+(b—a)W(t—h); etc. 


To compute the fundamental element of this solution, v(t), one must first find the roots x, of Eq. 
(52b). For the present case where p=5 the first 6 of these are: 


n: 1 2 3 4 5 6 
Xa: 0.5150; 1.2462; 2.0098 ; 2.7833; 3.563; 4.360. 


Aside from these roots, y(t) involves only the ordinary Bessel functions of zero and first order which 
are given in Hayashi’s tables. We set now explicitly: 


k/u=1.8; f=0.2; B=4.5X10- atmos.~'; 7r,,.=20 miles, b=0; c=3a; d=1.857a; e=1.429a. (61) 


The first set of constants gives: 
~ 


t= 1.6686  10~‘t (days); t;=0.1235; ty =0.1273; tz = 0.1387: t,=0.1490. (61a) 
For numerical plotting in terms of quantities of practical interest it must be finally observed that: 
Ay =yo0BAp; a=P,(O¥/OP) ro = uBQ/2rkh, (62) 


where / is the thickness of the sand. Choosing now h=130’ and Q=1.05 10° bbl./day (350,000 
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Fic. 4. The establishment of the “reservoir pressure’ on 
shutting in a well. 4p/(4p);=(pressure differential after 
shutting in) /(flowing pressure differential). 7=kt, uffr2 is 
the dimensionless time. 


Fic. 5. Theoretical pressure decline for the ‘‘East Texas” 
oil field, assuming it to be producing by a “water drive.”’ 
Crosses are observed average pressures for the west side 
of the field, which correspond to the “effective field pres- 
sures” or the ‘“‘well’’ pressures of the theory. The dis- 
crepancy at the age of about 670 days (Jan., 1933) is due 
to a field shut down in the last half of Dec., 1932, which 
was not taken care of in the calculations. 


bbl./day for the 120° sector subtended by the ‘East Texas’ field)’ it is finally found that Eqs. (60) 
may be rewritten as: 


Pe— Pw (Ib.) = 633.7 |log p+y/(t)} ete. 


The calculations are plotted in Fig. 5. The ordinates are the values of p,, corresponding to the 
“effective field pressure”’ as the well in the present problem is really supposed to be the equivalent of 
the whole field. The external pressure p, has been taken as the initial reservoir pressure of 1620 Ib. 
The abscissae represent the age of the field. 

As this practical problem has been presented only for the purpose of illustrating the general theory, 
it will be discussed no further except to remark that the observed performance of the ‘‘East Texas’’ 
field, up to September, 1933, follows, in a general way, the performance implied by Fig. 5 as may be 
seen from the crosses which give the observed effective pressures. The discrepancies can be largely 
accounted for by the fact that the production rates assumed for the periods in Eq. (57) were averages 
rather than exact values. The only arbitrary element in the analysis has been the compressibility, 8, 
which has been chosen so that at April 1, 1933, the observed and predicted pressures are in approxi- 
mate agreement. The fact that this value of 8 has a value of some 10 times that of ordinary water 
would imply that there is an equivalent of 43 percent free gas dispersed in the water horizon, which is 
in itself of some geological significance, although we cannot enter further into it here. 


(2) A single well in a closed reservoir" 


For this problem the results will be given only for the conditions of most practical interest which 
are, namely, that the flux at the external boundary r=,, is zero—a closed reservoir—and that the 
internal boundary radius is kept at a fixed pressure (derfsity). As Eq. (56) is valid, for any problem in 
which a flux is specified on one boundary while the density on the other, it may be adapted to the 
present case by setting: 


1° The author is indebted to Mr. R. D. Wyckoff for the '! The special case where f; is const. is the case of ‘“‘con- 
data used here for r., re, hand k/u which he derived froma _ stant well pressure” treated by Moore, Schilthuis and 
detailed analysis of well logs and of the flow mechanism ir Hurst, reference 9. 
the water reservoir adjacent to the ‘“‘East Texas’’ field. 
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wn 


re=r, (reservoir radius); f(t) =0, 
r1=?w (well radius) ; fill) =e, (63) 
g(r) =yi=const. 
Eq. (56) then takes for form: 


J (@ntw) J (ante) U(a,rye*2"* | 
Y=Vu -:- 0) ere a ty (64) 
JP (nh uw) —J (ant) 


for the density distribution within the sand at any time ¢. 


Of particular interest are the density decline at the outer reservoir boundary y(r,) and the decline 
in the flux from the well. In the notation of Eq. (42) these will be given by: 


(V¥i-— Yer) J (xn) Ji (xnp)e rntt 
Ye= Fot+2—__—____—_——_ (65) 
p Xn Je(xXn) —JSi?(xnp)} 
and 


4rkh(yi-—Yw) JP(xnpye xn2t 








: (a, (66) 
uB JP (Xn) — J? (Xnp) 
where x, is to be determined as a root of: 
Vi(xnp)Jo(Xn) — J \(xnp) Vo(X2 )= 0, (67) 


p having a value of about 2000. This large value of p makes the determination of the roots rather 
laborious, in view of the limited range of the available tables of Bessel functions. The numerical 
discussion of the above equations will, therefore, be omitted. 


(3) Well of infinitesimal radius 


When the flux at the well is specified one may replace the well by a mathematical sink and thereby 
simplify the analysis. Restating the problem in the form: 


lim (7 dy or) = Jos) 
r—>0 

v=fe(t) > %=TLe, | - (68) 
y=a(r) :t=0, 


the Bessel functions of the second kind must be excluded and the general solution is easily shown to be 
given by: 


2 Jo(anrye Kant 


t t 
7= yf re()Jolant)dr—« f fu(d)erer?™ dy + cant Ii(aat,) { firyenanl (69) 
‘ie JP (ante) 0 0 0 


where: 





J (a n%e) =()., (69a) 


For the special case where: f,,(¢) =const. = qo; f.(¢) =const. = y;; g(r) =const. = y;, Eq. (69) reduces 

to: . 
r 240 J o(anr)em —* 

Y= vito log —+ : 
re te ant T (ant) 








(69b) 


Although the roots of Eq. (69a) may be readily found from tabulations, the present case is not of 
very great practical interest and hence will not be treated here numerically. 
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D. Solution for Case (3) 
Again the functions (7;, y2, ys) will be introduced where: 


ria Or) Yi» ¥2> 73) =(1, 0, Q) -r=7)\, 


=(Q, 1,0) :r=r,; 


2, 


70) 


~ - A“ 


Yils Y2s 73 


= (0, 0, g(r)) : t=0, 


the complete solution being given in terms of (7:1, v2, ys) by Eqs. (32) or (51). For the Bessel function 
expansions the same form for U’(a,r) will be taken as in Eq. (52), but it will be required that the a, 


are the roots of: 
U"(a,r1) =0 (71) 


rather than of Eq. (52b). By using this definition and applying Eqs. -(27)—(30) once more to find the 
expansion coefficients, it is found that: 





(r?+4xt) re" 1 prmtr | 
ee ee re+3r.*)+r27(ri? log r,—r2? log re) 
Dogon) wort Gaol a 
Tv J (aunt1)J;(anre) Ula,rye Kant 
si Nae a ere —-, (72) 
ry an} J(anr1) —J7(a,re)} 
Peas) ont 1 pntor i | 
eid ae ; Pat Molen ro” + 3r;*) +—(r2? log re—r;" log ny 
2 ro — 1") i the A r2—r,2)?| 4 
° = P 
7 JP(a,r)L (a,r)e—xan't 
+—) — eee —_——_, (73) 
re Qn} J P7(a,7%;) —J;7(@,72)} 
P. : i a7? J;7(anr1) Ula,rye can*t T? 
re nd ot rg(r)dr+ —— ee ee a rg(r)U(a,r)dr 74) 
(fo°—7r;")e > ? J;? ant:)—J(a,re) r, 


so that finally: 
2k : 2 
y¥=—— J i fel(A)—fil(d) {dvA+ = rg(r)dr 
"y 


2 2 2 2 
re" — 1" 9 ts gimed dw. 


an? J (ant) U(a,re" Tr re 
+ 7) ——_-+—________ . --- ~J \(a@,7;) rg(r)U(a,r)dr 


JY (ant1) —J17(a,,1r2) 2 vr, 


J; ane) 7 . J (a7) v! a 
+ «—— f fi(dr)ese dd — x— { Fo(d)e*e" dd]. (75) 
TiQy, n Fo p 


Eq. (75) might be applied to the problem of the pressure decline in a field producing artesian water 
or an oil field driven by a water drive in which the. flux into the water reservoir, such as that due to 
rain fall, is known or assumed. Of more practical interest, however, is the case where the sand is 
closed off at the external boundary and the internal boundary represents an actual well of radius 


'? The constant terms in Eqs. (72)—(81) enter for case (3) aries. The: phenomenon is similar to the entry of the con- 
because of the vanishing of U’(ar) on both boundaries. stant term in the Fourier cosine series expansions. In this 
For Eqs. (77)-(81), the expansions are Dini series expan- connection it may be observed that the solution given by 
sions which take on a constant term for the zero order Carslaw (p. 118) for the special case of Eq. (81) in which 
Bessel functions when the derivative at the limiting radius. f-=f,=0, does not contain the constant term and will, 
vanishes. Eqs. (72)-(75) are generalizations of the Dini therefore, give the erroneous result that an insulated disk 
series expansions and take on constant terms when the _ will always tend to a zero temperature equilibrium regard- 
derivatives of the Bessel functions vanish at both bound- less of its initial nonvanishing temperature distribution. 
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small as compared to that of the external boundary. Here again, the well may be replaced by a 
mathematical sink, and for the more general case where the external boundary flux is not taken to be 
zero, i.e., Where: 


lim (7 dy/0r) = f(t); (r Oy Or) re = f(t); y=2(r) :t=0; (76) 


r—i() 


the functions (y1, ye, y3) take the form: 


3. (r®?+4xt) r 2_ JIflasrle 
"n=-- +g —+—7Z, (77) 
4 2r.* t, Pel and la, iP. 
—1 r+4«t 2 Jolanr)e~*%"! 
Y¥2=— + om 8 , (78) 
4 2r.? Pe nS (Ante) 
2 2 SIolan,rjen*a""*§ pv" 
¥3= { rg(r)dr+--——-), —_—_——_——_ f rg(r)Jol(anr)dr, (79) 
r*J Ps J (a,%) 0 
where 
J (apr) =0, (SO) 
so that 
2k vf 2 re 2 J y(anrye Kan? t re 
y= q] 1 fe) — fu A) {dX+ { rg(r)dr+—)>_ —_—— lf rg(r)Jo(a,r)dr 
Pe 0 redo Pu J (ant) 0 


t t 
—k { fa (dere r°AN + KS (Ante) { fuinyerran| (81) 


If the pressure, and hence density, is initially uniform (g(r) = y;) and the sand is closed off (f. =0) 
Eq. (81) takes the form: 


2p Tylanresantt at 
Y=Y7i- if fe(d) d+ —_——_—— f fu(nyeran], (81a) 
se 0 JP (ant) 0 


If now the flux at the well, f..(¢), has the constant value gq, y is given by: 


- J o(X ren 2"! 
v =vitd +log 7— 3(7°+4t) +2);-——_— | (81b) 
KS (x n) 
where: 
r=r/r: Xn=Qale: t=xt/r2. (81c) 





Even without any numerical assumptions it is that it will not affect the density, and hence 
clear that Eq. (81b) involves two distinct types pressure distribution. That is, the density 
of transient. The first, represented by the series, gradient will remain fixed, except for the series 
will disappear at an exponential rate as the time transient, although the absolute value of y will 
increases. The second, however, represented by decrease linearly with the time. 
the term linear in ¢ will persist indefinitely (since From a numerical point of view, however, the 
q has been assumed to be maintained indefinitely) series transient turns out to be of such short 
and gives a constant rate of variation of the 
density with the time, after the first transient 
has’ become of negligible magnitude. Further- 
more, this linear term is independent of r so k/u=1; f=0.2; B=4.5X10-4 atmos.'; r,=660' 


duration, for cases of practical interest, as to be 
of no importance. Thus assuming that: 
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Fic. 6. Pressure decline due to a single well at the center 
of a closed reservoir producing 1 bbl./day ft. of sand. 
Curve I: Pressure drop (from initial value) at well (1/4’ 
radius). Curve II: Pressure drop at external boundary 
(660’ radius); constants are those of Eq. (82). 


so that: _ 

t= 2.372t (days), (82) 
and applying Eq. (81b) one finds for the drop in 
pressure at the well (r7,.=}’) from its initial 
value, for a production rate of one bbl. ‘day per 
foot of sand, the expression: 


Pi- Pw 
- é i 
=0.1412( 7.12854 27-25 pineanio —), (83) 


x2 Se (x2) 


where Jo(x,7,.) has been replaced by 1, due to 
the small value of 7, (1/2640), p;— p. having the 
units of lbs. 

For the pressure decline at the external 
boundary of the sand, r=r,(7=1), it is found 
that: 


Pi- pe __ gs 
a e Tn*t 
=0.1412( -0.25-+27-25 — ). (84) 
XJ o(Xn) 


Eqs. (83) and (84) are plotted in Fig. 6. 
Because of the extremely short life of the series 
transient the time scale has been taken as hours. 
For as appears from the curves, after a period of 
only two hours the system settles down to an 
effectively steady state drainage by the well, and 
hence a uniform pressure decline over the whole 
reservoir, always maintaining a pressure differ- 
ential of 1.04 lb. between the reservoir boundary 
and the well. The very rapid pressure decline at 


























Fic. 7. Pressure distribution in closed reservoir as in 
Fig. 6 after the passing of the initial transient. r,=external 
boundary radius (660’); p= well pressure. 


the well and the lag in the decline at the external 
boundary are, of course, to be expected. 

The pressure distribution after the initial 
transient is clearly not a strictly steady state 
distribution, although the logarithmic term, 
characteristic of the strict steady state, does 
predominate everywhere except at the external 
boundary. Using the same constants as in Eq. 
(82), the pressure variation after the series 
transient has disappeared may be expressed as: 


p(r) — p. =0.1412(7.8785 —r*  2r,2—log r,/r). (85) 


Eq. (85) is plotted in Fig. 7. In general ap- 
pearance it resembles the strictly steady state 
logarithmic distribution and, in fact, approxi- 
mates it very closely as the well is approached. 
At the external boundary, however, the gradients 
vanish in Fig. 7, since the reservoir is closed, 
whereas in the steady state system they are of the 
order of 1/r,. As the result of these smaller 
gradients the total pressure drop p, — py = 1.04 Ib. 
is less than the corresponding value of 1.11 Ib. 
required to give the same production in an 
incompressible steady state system. 

One may nevertheless derive from these 
considerations a rather simple physical approxi- 
mation to the description of the decline in the 
system. For as has been seen above, after the 
series transient has passed, the density decreases 
uniformly over the whole reservoir and at a rate 
proportional to g. And, in fact, if Q be the actual 
production rate it is readily seen from Eq. 
(81b) that: 


O= frr2(d/at) = frr2yB(dp/at). 


(86) 
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Thus the production from the reservoir is sup- 
plied by the equivalent density or pressure 
decline distributed uniformly over the reservoir, 
the instantaneous pressure distribution differing 
only slightly from the steady state type. We are 
hence led quite naturally to the representation 
of the system as a “continuous succession of 
steady states’ already taken as a hypothesis in 
the theory of gas flow,* where the rigorous 
solutions could not be obtained explicitly. 


E. Non-radial flow; Green’s function 


Thus far we have treated problems involving 
single wells placed at the centers of their 
reservoirs. In practical cases, however, the 
reservoir will be covered by a number of wells, 
and the question arises as to the extent of their 
mutual interaction and their interference re- 
lations. Although we shall not enter here into all 
the details of the problem, we shall present the 
fundamental element by means of which the 
general solutions may be synthesized. This 
element is essentially the Green’s function of the 
circle which we shall define here as the function 
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y having a unit logarithmic singularity at (r’, 6’) 
with a vanishing normal derivative at the circle 
r=r, and which is identically zero at ¢=0. 
Physically, this Green’s function gives the 
density distribution within a closed circular 
reservoir, initially at ‘“‘zero density,’’ containing a 
well at (r’, 0’), i.e., displaced from the center of 
the reservoir. 

With the notation of Fig. 8, the constant 
strength well or sink may be represented by the 
term: 


vyo=log r’’/r.=3 log {r?/r2—2r'r/r2 cos 0+7r”/r2}. (87) 


To get a solution with a vanishing normal derivative (at r=r.) we add the function: 


¥a= 4 log (1—(27’r/r.2) cos 04+7r7r?/r.4) — (1/272) (7? +4 xt). (88) 


For our final solution, or Green’s function, we may, therefore, write: 


Y¥=VotvVatd A nmJn(Gnm’) Cos nOe~*@nm*t, (89) 
where 
J n/(tnm?e) =0 (90) 
and - 
vot Ya(t=0) +20 A nmJn(Qnmr’) Cos nO =0, (90a) 


nm 


so that Eq. (90) determines the parameters a,» and Eq. (90a) the coefficients A nm. Expanding yo 
and y, as Fourier series and equating first the coefficients of cos n@ we get: 


D A nmI n(Qnm?) = (r'r/r2)"/n+by 


m 


bn,=(r/r’)"/n :r<r'sn>0; (91) 


=(r’/r)"/n :r>r’ 


7 A omJ o( om?) => r’?/2r2+bo 


m™m 


bo =log r,./r’ i: r<r’ (91a) 


=log r./r:r>r’. 
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In virtue of the Bessel equation, recurrence relations, and Eq. (90) which the J,,(a@nnr) satisfy, it 
follows that: 


f tJ 1?(cnmr dr = (nmr? — n?) / Zon? JJ n?2(nmle)s (92) 
0 
i] ” 1D. (Qnmr)dr = (nr, _ Gan’) J a(Gan¥e)> (92a) 
0 
f{ b,rT n(Qnmr)dr —_ (1 Ginn?) 2In(A@nmt’) me, (r’/re)"J (Qnml ) }, (92b) 
0 ° 
f 43 J o(aomr)dr = (397 Gon?) J o(Qomle), (92c) 
0 
f bor J o(aomr)dr = [ Jo(aomr’) — J (om? ) | ‘om. (92d) 
0 
Setting now: 
Ganle™ Semi kt/r2=t; r/te =f; r'/r. =p, (93) 


Eq. (89) may be finally expressed as: 
y= +} log (7 — 27 cos (6’— 6) +?) +3(3 — p*?) +3 log |1—27p cos (6’—8) + p*r*} —3(7 +41) 


J o(X0mp) J o(XomF em 20m7* ad J (x ant St Afaat) ™ 
+2>>m +4>-n> m- cos n(0’—@)e~zm**#, = (94) 
Xom? J o?(Xom) 1 (x on” _ n*) J (x —_ 








This expression gives the density distribution in a closed unit radius (in the units of Eq. (93)) 
reservoir initially at density 0, which is drained by a unit strength well at (p, 6’). As we might 
expect, this Green’s function is symmetrical in (7, 6) and (p, 6’); that is, the density at (7, 0) ata 
time #, due to a source at (p, 6’) is the same as that at the same time at (p, 6’) due to a source at (7, @). 

If the initial density has a distribution given by g(r, @) the function: 


1 Qn 1 J o(Xom¥)e~z0m?= 1 Qe 
Yi=- { 7 d7 { g(7, 0)d0+-m { 7 I o(Xom?)d? f g(7, 6)d0 
Tre 9 J 5 us J 0?(Xom) =o 9 





2 Sau* J ASant) cos nOe-Znn*t 1 aa 
—_ f FIaCcmnt)dr f e(7, 8) cos 2n6d@ (95) 
0 + 


7w "70 (x an* —_ n*) J 1?(Xnm) 





should be added: 


If, finally, we suppose that the well has a flux strength g(¢) instead of the permanent unit value, 
the resultant density distribution will be: 


v=1F,0+ f g(d)(d/dt)yx(F, t—d)dd, (96) 


where 7; is the function of Eq. (94). 

By adding together a number of functions as in Eq. (96) with different values of r’, 6’ and q(t), 
the density and, therefore, pressure distribution due to the corresponding group or network of wells 
in a closed reservoir will be obtained. A study of such distributions will give the interference effects 
and mutual interactions among the wells due to their simultaneous drainage of the same reservoir. 
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Fic. 10. Pressure (density) decline at a central well as 


' é affected by the addition at ¢=0.3, of an off-center well 
Fic. 9. Diagrammatic representation of a two well inter- producing at same rate. (yi—7(1))/ ‘q=density drop at 


ference system. central well per unit flux. = ‘dimensionless time.’ 


As an example of such an interference problem and to indicate in more detail the method of 
treatment we shall consider the simplest case of interference, that of two wells in a closed circular 
reservoir. As is indicated in Fig. 9, we will consider well No. 1 placed at the center and well No. 2 
at a radius 7’, this radius being taken as the polar axis. We shall also suppose that No. 1 has been 
producing at a uniform rate of q; units since the initial instant and that No. 2 is opened after a time 


t; and is produced thereafter at the rate of g2 units, the production of well No. 1 being maintained 
at q1. 

Up to the time ¢; the pressure at well No. 1 will evidently decline without any interference, and 
according to curve I of Fig. 6 which corresponds to this case. The opening of No. 2, however, will 
clearly accelerate this decline, and it is this reaction of well No. 2 on No. 1 that we shall compute. 
One could also compute the interaction of well No. 2 on itself, i.e., the effect on the normal decline 
at the position of No. 2, as caused by the production of the central well, due to the drilling of No. 2; 
but to show the method the reaction of No. 2 on No. 1 will suffice. 

As just mentioned, for t<¢, the pressure at No. 1 (at r=r,,) will decline as in a radial system. 
Hence by Eq. (81b), to which Eq. (94) reduces for this case, we have for the density at (1): 

t=h 
¥(1) =yvitqil$—log r./re—2t+2E [em /x 2 I e2(xn) I), (97) 


where (7,72)? has been dropped and Jo(x,7./r-) has been replaced by unity, y; being the initial 
uniform density. 


For i? this contribution will continue but upon it will now be superposed the effect due to 
well No. 2, which by Eq. (94) is: 


go f3 Jo(Xomp)e~ 7” Pmt) 
y2(1) =qe log + (-— *)- 2q2 (t—t1) )+2q2 >—_——— a, (98) 


Xom” 5 Jen* xs.) 


This is to be added to the y(1) of Eq. (97). For the case where g:=qg2=4q, p= 3, and r,./r~. = 2640, we 
find for the density drop at well No. 1 the numerical expressions: 


t=t 
Cvi— (1) J/q=7.129+2t— 23 [en 27/2 IP (xn) J (99) 
th 


Cyi— (1) ]/q = 7.197 +2(28—fy) — 2D [e207 x2 Te2(xn) 1A To(xn/2)e7m}, (100) 
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These equations are plotted in Fig. 10 with ¢; chosen as 0.3. The dotted curve indicates the 
pressure decline if No. 2 well had not been drilled. It will be seen that after the initial transient, in 
which the original decline is altered, the rate of decline assumes a value twice as large as its original 
or extrapolated value. This is, of course, to be expected since for f>0.3, twice as much fluid is 
being taken from the reservoir as before. 


IV. Non-STEApDyY STATES; GAs FLow 


For the non-steady state flow of gases we must return to Eq. (7). As has already been mentioned, 
this equation is nonlinear and fo general solution for it has yet been found. We shall, therefore, 
outline an approximation method, which is similar to one applied some time ago to the theory of 
isothermal gas flow.* This is based on the physical assumption that the flow may be described by a 
continuous succession of steady states with the time entering as a parameter in the boundary con- 
ditions. The initial transient, however, which was neglected in the earlier solution to the isothermal 
gas flow problem will be taken care of by introducing a “‘radius of drainage’ which recedes to the 
boundary of the reservoir as the steady state condition is established. That is, the “radius of drain- 
age”’ gives the distance from the well at which the approximately steady state distribution has been 
established, and the production during the initial transient is simply the rate at which gas is removed 
from the system in the establishment of this steady state distribution, it being assumed that no 
gas is removed from any point until the radius of drainage has passed that point. 


mA 


Fic. 11. Diagrammatic representation of the initial transient in a gas flow system. ro=“‘radius of drainage.”’ For r>r», 
system is undisturbed; for r<rpo steady state distribution is already established. 


Considering now the specific problem of a closed reservoir of radius 7, and a central well of radius 
r., we apply the steady state pressure distribution of Eq. (14a). Thus assuming that the initial 
reservoir pressure is p; and that the well pressure is permanently maintained at p,., we will have for 
the production rate, at a time ¢ when the “radius of drainage” has receded to the radius ro (the 
period of the initial transient) the value (cf. Eqs. (14), (15); and Fig. 11): 


2rkAy"t™!™ rf) re re 
Q =—____—_——_——__ == -2 rf Lf vrdr+ f virdr| (101) 
(1+) pyo''” log ro 7, OtlLe ry r 


vyo=7(p=1); Ay “t= my (itm ma oy, (Ite) / me y¥i=y(t=0). (101a) 





where 


By introducing the notation: 


2 log 70/fu =2; 2 log r/r.=¥, and applying Eq. (14), (101b) 








Is 
al 
or 


a) 
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Eq. (101) may be rewritten as: 





myo "ufr,.? pidz pr? yerdy 
Pet. f om (102) 
0 


4k g Jo [(Ayot™ m 2)yty¥et™ all (+a) 


In principle, Eq. (102) permits the determination of ¢ as a function of z, or of ro as a function of ¢ 
and hence the time variation of Q, on application of Eq. (101), till ro recedes to the reservoir radius 
r,. To show the nature of the process it will be sufficient to choose the case where p.=0 (=) 
when Eq. (102) can be reduced, by partial integration to the form. 


t=[m(1+m)u(yo/ve) "fru sei fev—eutul 2)" i dy, (103) 


which is rather convenient for graphical integration. 
Taking the special case of isothermal flow, m=1, Eq. (103) gives, when evaluated graphically: 
ce? =47o/r7~=1+1.068(4Rp,t/ fr,,2)% 9 =2 /r,.(Rpt/ pf) (104) 
except when /~0. 
With the constants: 


f=0.2; k=1 darcy; u=0.012 centipo'se; p. = 100 atmos; r,o=1/4’. (104a) 


Eq. (104) is plotted in Fig. 12 as curve I. It will be observed that the ‘‘radius of drainage’’ recedes 
extremely rapidly, growing to a radius of 100’ in about 1 minute and reaching the effective reservoir 
boundary (r,~660’) in less than an hour. 

The production rate which is given during this initial transient by: 

Q=2rkyop.?/uz = 3518/2 g/sec./cm = 4.86 X 108/z ft.*/day/ft. (sand) (105) 
is plotted as curve II with the constants as in Eq. (104a) and yo taken as 6.72X10-* g/cc, the 
atmospheric density of methane. 

After the initial transient has passed, i.e., after the ‘‘radius of drainage” has reached the reservoir 


boundary we shall assume that reservoir pressure will decline uniformly throughout, its instantaneous 
values following a steady state distribution. Thus we shall now have instead of Eq. (101): 








4mrkAy tm! 0?” — B70 fy.''™ Oye 7°" yerdy 
Q =-——_——_ = -—-2rf- yrdr = -———_—_ --- ] ———- (106) 
(1 +m) psyo'!™ Ot rw z at Jo [Ay tm my ‘s+ V0 Vm) lem im 
with the notation of Eq. (101b). 
Again choosing the case, y,,=0= p,, we find: 
OY. —Aks 1/(1+m) 
% Ce ee en P (106a) 
ot . 
(1 +m)ufr.?yo' mf ym tmeudy 
Ai) 
=—2k/fuyosr?; m=1; r,/Vw = 2640, (106b) 
so that, with the constants of Eq. (104a), m=1, and recalling Eq. (5) 
p.=1/[1.129 X 10-2(t—to) +1/ pi], (107) 


where ¢=¢o is the time (in days) when the radius of drainage reaches the reservoir boundary. 
The production rate for the present case (m=1, p,.=0) is given by: 


O= thyop2/(u log re/fw) =2.23 X10 p? g/sec./cm = 3.09 X 10°p,? ft.3/day/ft. sand. (108) 
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Fic. 12, Production rates and recession of ‘‘radius of 
drainage’ during the initial transient about an isothermal 
flow of gas well. Curve I: Extent of radius of drainage in ft. 
Curve II: Production rate in 10° ft.*/day. Constants as in 
Eq. (104a). 
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Fic. 13. Reservoir pressure and production rate decline 
after the initial transient in a closed gas sand producing 
through a central well under isothermal conditions. Curve 
I: Pressure at boundary of reservoir (atmos.). Curve II: 
Production rate in 2.5 10° ft.8/day per ft. of sand. Con- 
stants as in Eq. (104a). 


Eqs. (107) and (108) are plotted in Fig. 13. Here again, the decline is very rapid, the pressure 


dropping to a tenth of its initial value in about 83 days." As to the relation of the above decline 
curves to those as are observed in practice, it should be remembered that not only is the assumed 
permeability of one darcy higher than that of most gas sands, but of more importance, the effective 
reservoirs drained by each well in a field will be considerably higher than that of 660’ assumed here. 
Both these factors will tend to increase the lives of the wells over that indicated by Fig. 13. 

Finally, it may be observed that the above theory may be readily extended to cases where the 
well pressures are not zero, and where the flow is not isothermal. In principle such cases are all 
covered by Eq. (106), although the numerical integrations become more complex. 

The author is indebted to Mr. R. D. Wyckoff for many helpful discussions of the problems treated 
in this paper and to Dr. P. D. Foote for permission to publish it. 


‘8 The apparent discrepancy between this high decline non-steady state gas flow, to that given in the earlier paper, 


rate and the much slower one previously derived by 
Muskat and Botset’ is due to a different choice in the 
physical and geometrical constants. For although the 


from a physical point of view, the numerical conclusions for 
the two theories are almost identical, for the same con- 
stants. 


author prefers the present approximation theory of the 
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Sorption on Cotton Fibers of Dyes with Varying Molecular Association in Solution 
GEORGE L, CLARK AND JULIA SOUTHARD, Department of Chemistry, University of Illinois 
(Received December 9, 1933) 


X-ray diffraction patterns of aqueous solutions of the 
oxazine dye, Nile blue sulfate, yield direct evidence of a 
varying molecular association as a function of concen- 
tration, which was previously indicated by Cohen's 
potentiometric and spectrophotometric data. Solutions of 
the order of millionth molar approach obedience to the 
laws of dilute solutions; in moderate concentrations of the 
order of 5X10~ molar the dye molecules associate to 
micellar structures; and in concentrated solutions the 
association proceeds to the stage of precipitation. The 
largest change in the value of the identity spacing d 
corresponding to the inner edge of the principal diffraction 


I. INTRODUCTION 


HE problem of the mechanism of sorption of 

dyes on textile fibers still retains great 
interest, not only from the practical point of view 
of successful commercial dyeing of fabrics, but 
also because of the possible light thrown upon 
the fine structure of a given fiber by its behavior 
in the dyeing process. A pertinent example is to 
be found in surgical absorbent cotton. What 
structural factors are involved in producing 
cotton with maximum absorbency is a question of 
immediate significance. Obviously the sorption 
of dyes on any fiber depends on the basic or 
acidic nature of the dye and on the polarity of the 
surface of the absorbent.! The porosity and the 
size of the capillary pores in the fiber must also be 
involved. It has been demonstrated? that water is 
soaked on cellulose in two stages. First it 
attaches itself to molecules of cellulose; and 
second, it enters the capillary spaces. The theory 
is that water or other polar molecules such as 
dyes attach themselves to the alcohol groups, one 
primary and two secondary, in cellulose. By use 
of a formula relating vapor pressure and curva- 
tive of surfaces Urquhart and Williams’ calcu- 


'S. E. Sheppard, R. H. Lambert and R. L. Keenan, J. 
Phys. Chem. 36, 174 (1932). 

*S. E. Sheppard and P. T. Newsom, J. Phys. Chem. 34, 
1158 (1930); 36, 930 (1932). 

*Urquhart and Williams, The Moisture Relations of 
Cotton, Shirley Institute Memoirs III, 197 (1924). 
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halo is found to begin with about 10~* molar solutions. The 
attempt is then made to measure pore sizes in cotton fibers 
from the relative absorbability of the dye from its solutions 
over such a range of molecular and aggregate dimensions. 
Nile blue sulfate and methylene blue solutions before and 
after addition of absorbent cotton were subjected to 
spectrophotometric analysis. Again it is indicated in the 
case of the former that there is a large increase in dye 
absorption by the cotton corresponding to the change with 
increasing dilution from micellar to molecular dimensions 
which are compatible with predominating cotton pore sizes. 


lated the average diameter of the smallest pores 
of cotton to be 5.110-* cm, or about 1} times 
the diameter of the water molecule. These 
measurements indicate that only relatively small 
molecules can really enter the smaller pores in 
cotton fibers. 

It is conceivable that a whole series of dyes, 
for example, representing a progression in 
molecular sizes, might be found with which pore 
structure might be ascertained by sorption. 
However, necessary chemical differences in the 
series would introduce a variable such as to 
constitute an insuperable obstacle in the way of 
dependable interpretation. One might despair of 
realizing the ideal sorption condition, namely a 
single substance whose molecular size in solution 
could be varied as a function, for example, of 
dilution. Such a dye would make it possible to 
determine at once the range of concentration and 
of molecular or aggregate size over which maxi- 
mum sorption on cotton took place, independent 
of the variables of chemical constitution and 
characteristic specific absorption. These data 
would then give an approximate measure of 
compatible pore sizes. Singularly; enough a dye 
does actually exist which fulfills these ideal 
requirements. It is the purpose of this paper, 
therefore, to present the results of experiments 
with this dye, both as to the effects of concen- 
tration in aqueous solution upon molecular 
aggregation and the sorption on cotton fibers. 
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Il. Tne MoLtecuLaR AGGREGATION OF NILE 


BLUE SULFATE IN SOLUTION 


The remarkable properties of the oxazine dye, 
Nile blue sulfate which is widely used in his- 
tological and biochemical practice were first 
demonstrated by Cohen and Preisler* in the 
course of their beautiful series of studies on 
oxidation-reduction. Both from e.m.f. and from 
spectrophotometric data, it was demonstrated 
that there is a very marked concentration effect. 
For example, there is an absorption peak fixed at 
630 my while another changes from 578 at a 
concentration of 0.0018.17 to 595 at 0.000009.1/, a 
shift which is approximately a function of the 
logarithm of the total concentration. Hence the 
molecular species or tautomer predominant in 
concentrated solution changes its structure in 
step wise fashion with dilution until finally it 
becomes identical with the presumably simple 
structure characteristic of dilute solutions. The 
only type of tautomer that can answer this 
description seems to be an associated molecule or 
a micelle whose complexity varies with the 
concentration of solution. The evidence is 
summed up as follows: 

(1) Nile blue in concentration of the order of 
millionth molar in water approaches obedience to 
the laws of dilute solutions. 

(2) In moderate concentrations (110~° to 
5 <10~*) Nile blue acts as if it associated to form 
micellar structures. This takes place in acid 
solutions down to pH 7. 

(3) In concentrated solutions the association of 
Nile blue molecules proceed to the stage of 
precipitation. 

As a test of these observations it seemed of 
interest to the writers to study Nile blue 
solutions by the x-ray diffraction method. 
Krishnamurti® was the first to show that the 
measurement of the intense ‘‘amorphous”’ scat- 
tering of x-rays at small angles by colloidal 
solutions leads to evaluations of aggregate or 
micellar size and that molecular weights calcu- 
lated data give results in 
agreement with other methods. Usually a diffused 


from diffraction 


‘Cohen and Preisler, Public Health Reports, Supplement 
92, 1931. 
* Krishnamurti, Ind. J. Phys. 5, 473 (1930). 


AND 





JULIA SOUTHARD 


broad halo characterizes the diffraction pattern. 
When the ultimate particles in solution have a 
size distribution, the inner edge of the diffraction 
halo is a measure of the size of the largest 
particles. 

The x-ray diffraction photographs were made 
with rays from a molybdenum target x-ray tube, 
in preference to copper radiation in order to 
produce as sharp haloes as possible. 

A small piece of sheet brass, 1/8 inch in diame- 
ter and 3/4 inch square, with a hole in the center 
covered with mica served to hold the solutions. 
Picein served to paste the mica windows to the 
cell. Mica was used as it gives no diffraction rings 
of its own but instead spots which are easily 
distinguished. A small slot was cut in one end of 
the cell to permit entrance of the solution. The 
cell was attached directly to the pinhole with the 
use of plasticene. 

Pure Nile blue sulfate was obtained from Dr. 
Barnett Cohen of the National Institute of 
Health. 

The following procedure was used in preparing 
the solutions: A sample of the Nile blue sulfate 
was carefully weighed out and dissolved in 100 cc 
of distilled water to give a molar concentration of 
0.0884. Ten cc of this solution was then taken and 
diluted to 100 cc, ete. Before each dilution 
was made, care was taken to be certain that all of 
the Nile blue sulfate was in solution. Solutions by 
this means were obtained 
concentration of 0.0884, 
0.0000884, and 0.00000884. 

A sample of each of these solutions was placed 
in a cell made as described above, through the 
slot, made for this purpose, was then covered 
with paraffin to prevent evaporation during the 
exposure. A pattern of water and one of the pure 
solid Nile blue sulfate were then taken in the 
same manner. The patterns so obtained were 
then carefully analyzed by means of micro- 
photometer curves. 


molar 
0.000884, 


having a 
0.00884, 


From a measurement of the microphotometer 
curves, the angle of diffraction for the inner 
margin, the point of maximum intensity, and the 
outer edge of the halo were found. 

A graph was then made plotting d, as found 
from the Bragg equation, for the inner edges of 
the halo against the log of the concentration. 
This is shown in Fig. 1. 

















EE mee 





SORPTION 





Change in Maximum 
Porticle Size with 
Chonge in Concentration 

















| Re 
| \ 








SN AU 





N 






































2949645 394645 $ 94645 594645 694645 


Log Concentration 


lic. 1. 


Taste 1. Nile blue sulfate solutions. Constants found from 
points of maximum intensities. 


d0.71 
Conc. log conc. 0 Sin @ din ALU, 
0.0884 2.94645 6° 16’ 0.109 3.26 
00884 3.94645 6° 8.5’ 107 3.31 
00884 $.94645 6° 19’ .110 3.22 
0000884 5.94645 ge BY 111 3.19 
OCCOO884 6.94645 6 27.9’ VM 3.19 


TABLE I}. Nile blue sulfate solutions. Constants found from 
radius of inner ring of the corona. 


0.71 
Conc. log cone. 0 Sin@ dinA.U. 
0.0884 2.94645 1° 34’ 0.0273 13.00 
00884 3.94645 1° 49’ 0317 11.19 
000884 1.94645 .0346 10.26 
OCCL884 5.94645 2 we .0436 8.14 
00000884 6.94645 Zz Se .0500 7.10 


TABLE III. Nile blue sulfate solutions. Constants found from 
radius of the outer limit of the corona. 


0.71 
Conc. log conc. 6 Sin @ din ALU. 
0.0884 2.94645 
00884 3.94645 12° 41’ 0.219 1.61 
000884 1.94645 12° 24’ 214 1.65 
6000884 9.94645 ia ¢ .209 1.69 
00000884 6.94645 i 197 1.80 
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All data are found in Tables I to II]. There are 
several points of unusual interest which appear 
from these data as follows: 

(1) The value of d corresponding to the 
maximum intensity of the solution haloes checks 
closely with the value of the principal liquid halo 
for pure water, namely 3.19A. This value is 
observed for the two most dilute solutions, 
increases to 3.22 for 0.000884.\/ and to 3.31 for 
0.00884.17 solutions and then decreases again to 
3.26 for the most concentrated solution. 

(2) A suspension of the dye is found for the 
0.088417 solution. The fact that this is a truly 
crystalline suspension is shown by the pattern in 
Fig. 2 which is markedly different from patterns 





Fic. 2. Typical diffraction pattern for concentrated solu- 
tion of Nile blue sulfate (0.0884 molar). (Spots due to mica 
used as windows in cell.) 


for more dilute solutions (Fig. 3) in showing two 
quite sharp rings. These same rings are the 
strongest interference maxima for solid Nile blue 
sulfate with spacings of 13.00 and 6.45A. The 
inner sharp ring seems to be related to the inner 
edge of the halo for the other solutions, which 
moves to smaller and smaller angles as the 
concentration of dye is increased, until between 
limits of 0.0884 and 0.00884 molar a great change 
occurs in the state of the dissolved Nile blue 
sulfate which is evidently changing from a true 
colloid to a crystalline suspension. 
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ic. 3. Typical diffraction pattern for dilute solution of 


Nile blue sulfate (0.0000884 molar). 


3) Another large change observed from the 
inner edge of the halo occurs between 0.000884 
and 0.0000884.1/ solutions as shown by the curve 
in Fig. 1. The solute evidently changes from the 
colloidal state and approaches the state of a true 
solution. This observation then is a 
of Cohen 


molecul i! 
and Preisler’s 


electrometric 


direct confirmation 


interpretation from their and 
spectrophotometric results. 

4) The outer edge of the corona probably has 
less significance but it is interesting to note that 
there is a progression towards smaller angles and 
larger values of d with increasing dilution due 
either to hydration in true solution of solute 


particles, or to some tautomeric shift. 


[I11. SoRPTION STUDIES OF NILE BLUE SULFATE 
\ND METHYLENE BLUE ON COTTON FIBERS 


It is easily apparent from the foregoing data 
that the degree of molecular aggregation in 
:queous solutions of Nile blue sulfate varies as a 
function of concentration. The relative amount 
of dye adsorbed from solutions of various concen- 
trations by cotton fibers should be, therefore, an 
indication of threshold capillary pore size in the 
cotton. Several methods might be employed for 
such a study but the simplest and perhaps most 


accurate appears to be a spectrophotometric 
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alter 


and 
sorption to the equilibrium state on known 
definite weights of cotton. 


measurement of solutions before 


For this work the cotton used was furnished by 
the courtesy of the Experiment Station, College 
Station, Texas. The fibers were removed from 
the seeds by hand. Fats, waxes and resins® were 
removed by 48 hours extraction with chloroform 
hot Soxhlet. The 
removed from the cotton by placing it in a 


in a chloroform was then 


vacuum desiccator for three days at room 


temperature. Samples for experiment were 
prepared from the extracted and dried cotton, 
Most of the shorter fibers had been removed by 
stapling. This gave samples composed of fibers of 
more uniform quality (growth and development). 
The cotton fibers which were to have Nile blue 
sulfate sorbed on them were dried in an oven at 
90°C, cooled in weighing bottles, and weighed 
ready for use. 

The Nile blue sulfate was weighed, made to a 
known volume, and the molality calculated as 
before. Part of this original solution was diluted 
ten times and another part 100 times in order to 
the 
molalities as those found to have characteristic 


give solutions of approximately same 
properties by Cohen and Preisler. Fig. 4 gives 
absorption curves for each of the solutions as they 
were made. A 25-ce portion of each of the three 
dye solutions was measured. Dye was sorbed 
from these solutions by the weighed cotton 
sample for one hour. Again curves were made for 
each solution. A Keuffel and Esser spectro- 
photometer was used. Readings were made on 
each dye solution every ten millimicrons from 440 
to 600 mu and every 20 millimicrons from 600 to 
700 mu. The percent transmission was plotted as 
ordinate against the wave-length as abscissa. The 
probable error for an average of five readings at 
any one wave-length between 480 and 700 my is 
+0.8 percent of the reading, while below 480 my, 
it is +1.2 percent. It was necessary to use cells of 
different depth for the transmission of light 
through solutions of such great concentration 
differences. 

Curves for the methylene blue were obtained 


*P. H. Clifford, Lucy Higginbotham and R. G. Fargher, 
Chemical Analysis of Cotton, Shirley Institute Memoirs III, 
31 (1924). 














SORPTION 





Vit E BLUE SULPHATE 


napiia Concentration #5x10°*m 
Thickness of So/jvtion /mm 
Concentration 45 .x/0%m 
Thickness of solution 85mm 
_——-- Concentration 1/2 x10 4m 
Thickness of Solution 509mm 
---——. Concentration #5 x 10°'m 
Thickness of solution /00 nM 








100 


90 





80 


70 


60 


50 


40 


30f- 


20 


Per cent Transmission 

















$20 
Wavelength in m fa 


560 600 


Fic. 4. Typical spectrophotometric curves for solutions 
of Nile blue sulfate. For each pair of curves, the lower and 
upper represent, respectively, the solutions before and 
after sorption by the same weights of cotton. 


by the same procedure as for Nile blue. The data 
for methylene blue are shown in Fig. 5. Fig. 6 was 
obtained by plotting difference in percent of 
absorption of light, before and after removal of 
dyed fibers at the point of maximum absorption, 
against the concentration of the dye. This 
indicates that a larger percent of dye is removed 
from the more dilute solutions. This is especially 
true of Nile blue. The evidences of molecular 
aggregation as a function of concentration are 
therefore borne out by the absorption of dye by 
cotton fibers. 

An x-ray pattern of cotton which had a large 
amount of dye sorbed upon it indicated that the 
dye has in no way altered the cellulose pattern, 
nor is there any evidence of the presence of the 
dye. 


DISCUSSION OF RESULTS 


1. The fibers which were introduced into the 
methylene blue bath were wetted much more 
uniformly and quickly than were those in the 
Nile blue sulfate bath. This behavior is attributed 
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Fic. 6. 


to the difference in properties and difference in 
structure of the molecules. 

2. Fig. 4 shows the results of the Nile blue 
sorption. The original solution, curve No. 1, 
shows a decided peak of maximum absorption of 
light at 580 my, and a slight peak at 640 mu. The 
solution No. 1’ was a portion of solution No. 1 
which had had dye sorbed from it. The curve 
shows a difference in maximum absorption at 580 
mu of 10 percent on the spectrophotometer. 
Solution No. 2 shows a point of maximum 
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absorption of light at 590 to 640 mu, while 
solution 2’ has a peak of maximum absorption at 
640 mu. There is a difference of maximum 
absorption at 620 mu of 39 percent with a shift of 
40 my toward the red. Solution No. 3 has two 
definite peaks, one at 590 and one at 640 mu. 
Solution No. 3’ shows three points of maximum 
absorption—550, 590, and 640 mu. The largest 
difference in percent of absorption was 49 percent 
at 640 mu. This very great difference in ab- 
sorption of light may be due to the fact that the 
dye solution of this concentration behaves as a 
true solution. The dye particles are probably of 
the order which can enter most easily the spaces 
of the cotton fiber. 

Curve 1, Fig. 6, shows that there is a decided 
increase in percent of absorption of Nile blue 
sulfate at about a concentration of 5x10~° 
mole /liter. Even when the dilution of the dye is 
increased many times, there remains a peak of 
maximum absorption at 590 and at 640 mu 
(630 given by Cohen). If the peak of maximum 
absorption represents different tautomeric forms 
or molecular aggregates of the dye, the cotton 
seems to sorb approximately equally from the 
equilibrium forms sorbing at 580 and at 640 mu 
with the balance in favor of the latter. 

Curve 2, Fig. 6, is for the methylene blue and 
shows a point of maximum change in sorption of 


the dye at 5X10~*. Methylene blue in concen- 


trated solutions shows two points of absorption 
600 and 660 mu. The point at 600 my gradually 
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disappears and the one at 660 my becomes more 
pronounced on dilution. There was a maximum 
change of 41 percent of absorption in the 5.1 
< 10-4 molal solution at 600. 

Here it seems likely that an equilibrium exists 
between tautomeric forms. The dye solution on 
dilution takes on a slightly violet tint. The color 
changes suggest a shift from one form to another. 
This dye may exist as either a quaternary 
ammonium salt or as a sulphonium salt. 


CONCLUSIONS 

1. The two basic dyes, Nile blue sulfate and 
methylene blue, each gave similar but charac- 
teristically different absorption spectra in the 
visible, and it was possible to compare the 
sorption of dyes from measured volumes on 
known quantities of cotton fibers. 

2. Diffraction patterns of the cotton made 
before and after sorption of Nile blue sulfate, 
gave no indications of any change in the cellulose 
pattern. 

3. Cotton fibers showed a decided increase in 
the percent sorption at the concentration where 
the dye particles showed a large decrease in size. 

4. It is highly probable that two or more forms 
of the dye molecule exist in aqueous solution. 
Absorption curves bear out this statement. The 
more concentrated solution of Nile blue sulfate 
represents the concentration at which molecules 
aggregate into particles of colloidal size, as 
demonstrated by x-ray diffraction data. 











